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DEVELOPMENT OF A MULTIPLE-PARA^-^ETER NONLINEAR 
PERTURBATION PROCEDURE FOR TRANSONIC 
TURBOMACHINERY FLOWS: PRELIMINARY 

APPLICATION TO DESIGN/OPTIMIZATION 
PROBLEMS 

Stephen S. Stahara, James P. Elliott, 
and John R. Spreiter 

Nielsen Engineering ^ Research, Inc. 
Mountain View, CA 


SUMMARY 


An investigation was conducted to continue the development 
of perturbation procedures and associated computational codes 
for rapidly determining approximations to nonlinear flow 
solutions, with the purpose of establishing a method for 
minimizing computational requirements associated with parametric 
design studies of transonic flows in turbomachines. The results 
reported here concern the extension of the previously-developed 
successful method for s ingle -parameter perturbations to simul- 
taneous multiple-parameter perturbations, and the preliminary 
application of that multiple -parameter procedure in combination 
with an optimization method to blade design/optimization 
problems . 

In order to provide as severe a test as possible of the 
method, attention is focused in particular on transonic flows 
which are highly supercritical. Flows past both isolated blades 
and compressor cascades, involving simultaneous changes in both 
flow and geometric parameters, are considered. Comparisons with 
the corresponding 'exact' nonlinear solutions display remarkable 
accuracy and range of validity, in direct correspondence with 
previous results for single-parameter perturbations. Initial 
applications of the perturbation method combined with an 
optimization procedure demonstrate the ability of the multiple- 
parameter method to work accurately in a design environment 
and establish its potential for reducing the computational work 
required in such applications by an order of magnitude. 
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1 . INTRODUCTION 


With the continuing success of advanced computational 
methods to determine solutions to increasingly complex fluid 
dynamic phenomena, it has become clearly apparent that 
in order to employ these methods in applications requiring 
routine high-frequency use, a means must be found to reduce the 
computational demands necessary in their straightforward 
application. While this need exists across a spectrum of 
aerodynamic users, it is particularly high in turbomachinery 
applications. There both the basic aerodynamic computation is 
time-consuming and the number of variable flow and geometry 
parameters are large, making any turbomachinery parametric or 
design study computationally expensive under the best of 
circumstances, and in many instances using more advanced codes, 
prohibitively so. 

The ultimate objective underlying this study is to develop 
the means of reducing substantially the overall computational 
requirements necessary for turbomachinery design or parametric 
studies by minimizing the actual number of "expensive" numerical 
flow solutions required. That such procedures are achievable 
has been successfully demonstrated in the previous phase (ref. 1) 
of this study. In that work, a perturbation method was developed 
and tested on a large number of nonlinear flow problems involv- 
ing single-parameter changes of a variety of flow and geometric 
parameters. Subcritical and supercritical flows past isolated 
blades and compressor cascades were considered, with particular 
emphasis placed on supercritical transonic flows which exhibited 
large surface shock movements over the parametric range studied. 
Comparisons of the perturbation predictions with the corresponding 
’exact' nonlinear solutions indicated a remarkable accuracy and 
range of validity of the perturbation method. 

The work reported here describes the continued extension 
and refinement of that perturbation technique, and has focused 
on the development of its capability for actual application 
to practical turbomachinery design problems requiring the highly 
repetitive use of computational codes to determine a large 
number of related flow solutions. Two primary tasks were 
involved. The first consisted of the extension of the method 
to treat simultaneous multiple -parameter perturbations. The 
second involved the combination of the perturbation method with 
an optimization procedure, and the preliminary application of 
that combination to blade design/optimization problems. The 
nature of the present work is both exploratory and developmental 
in that aspects of the procedure- - such as its validity, range 
of accuracy, and computational economy for multiple -parameter 
perturbation problems and its workability in an optimization 
design environment - -will be investigated, and a computational 
code for multiple-parameter perturbations will be developed. 
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2 . 


ANALYSIS 


2.1 Perturbation Concept and Previous Applications 


The classical approach of performing a perturbation analysis-- 
consisting of establishing and solving a series of linear 
perturbation equations - -appears an an obvious choice for the 
current applications. However, results from the initial phase 
of this study (ref. 2 ) demonstrate that for applications to 
sensitive flows such as occur in transonic situations, the 
basic linear variation assumption fundamental to the technique 
is sufficiently restrictive that the permissable range of 
parameter variation is so small to be of little practical use. 

An interesting alternative to the linear perturbation equation 
approach has recently been successfully examined in which a 
correction technique is used that employs two or more nonlinear 
base solutions. For that method, the basic perturbation 
solution is determined simply by differencing two nonlinear 
base flow solutions removed from one another by some nominal 
change of a particular flow or geometrical quantity. A unit 
perturbation solution is then obtained by dividing that result 
by the change in the perturbed quantity. Related solutions are 
determined by multiplying the unit perturbation by the desired 
parameter change and adding that result to the base flow 
solution. This simple procedure, however, only works directly 
for continuous flows for which the perturbation change does 
not alter the solution domain. For those perturbations which 
change the flow domain, coordinate stretching is necessary to 
ensure proper definition of the unit perturbation solution. 
Similarly, for discontinuous flows, coordinate straining is 
necessary to account for movement of discontinuities due to 
the perturbation. 

In a number of recent applications of the method (refs. 1-7), 
results have been obtained which demonstrate the accuracy, range 
of validity, and versatility of perturbation methods based on 
such ideas. The most extensive and systematic of these are 
provided in references 1 and 7, where results are reported for 
a variety of flow and geometry parameter perturbation case 
studies of nonlinear subsonic and transonic flows past both 
isolated blades and compressor cascades. In those results, 
particular emphasis was placed on strongly supercritical tran- 
sonic flows which exhibit large surface shock movements over the 
parameter range studied. Comparisons of the perturbation results 
with the corresponding ’exact' nonlinear solution display 
remarkable accuracy across the spectrum of examples studied. 

The basis of this accuracy lies in the use of coordinate 
straining. This provides the means in determining the unit 
perturbation to account properly for the displacement of 
discontinuities and maxima of high- gradient regions due to a 
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parameter change. This in turn enables the perturbation method 
to maintain uncommon accuracy in regions of high gradients over 
large parametric ranges. In reference 1, a detailed examination 
was made of the effect of employing different classes of straining 
functions. Those results have illustrated deficiencies in 
certain classes of straining functions and have lead to the 
identification of a superior class of straining functions. These 
results are discussed in more detail in section 2.3. 


2.2 Simultaneous Multiple-Parameter 
Perturbation Formulation 

To provide the theoretical basis of the perturbation 
method as applied to simultaneous multiple-parameter perturbations 
of flows containing multiple shocks or high- gradient regions, 
consider the formulation of the procedure at the full potential 
equation level, since all of the results presented here are 
based on that level. Denote the operator L acting on the full 
velocity potential $ as that which results in the two-dimen- 
sional full -potent ial equation for i.e.. 


L[<D] = 0 (1) 


If we now expand the potential in terms of zero- and higher- 
order components in order to account for the variation of M 
arbitrary geometrical or flow parameters q- 


M 

$=$ + +... 
O . J 1 . 

J=1 J 




( 2 ) 


and then ins 
the result, 
components , 

e . - Aq . we 
3 ^3 

zero- and M 


ert this into the governing equation (1) , expand 
order the equations into zero- and first-order 
and make the obvious choice of expansion parameters 
obtain the following governing equations for the 

first-order components 


L[$^] = 0 




(3) 
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Here is a linear operator whose coefficients depend on zero- 


order quantities and 9L[5>^]/9qj represents a 
to the perturbation. Actual forms of L- 


forcing ' 
and the 


term due 
forcing ’ 


term are provided in reference 2 for a variety of flow and 
geometry parameter perturbations of a two-dimensional turbo- 
machine, and in reference 4 for profile shape perturbations of 
an isolated airfoil. An important point regarding Equation C^) 
for the first-order perturbations is that these equations 


represent a unit perturbation independent of the actual value of 
the perturbation quantity z . . 


Appropriate account of the movement of a multiple number of 
discontinuities and maxima of high- gradient regions due to the 
perturbation is now accomplished by the introduction of strained 
coordinates (s,t) in the form 


where 


X = s + 
y = t + 

x^(s,t) = 
y^Cs^t) = 


M 

I £ .X. (s , t) 

j=i ^ 

M 

I e.y. (s,t) 
j=l ^ 


N 

I 6x^x^ Cs,t) 
i = l i 

N 

I Sy-y^ Cs,t) 
i = l i 


(4) 


(5} 


and e.6x., e . 6y . represent 

J ^ ^ th 

perturbation of the qj p 

and x^ (s,t), y^ Cs,t) are 


individual displacements due to 
arameter of the N strained points, 
straining functions associated with 


1 1 

each of the N strained points. Introducing the strained 
coordinate Equations (4) and (5) into the expansion formulation 
leaves the zero-order result in Equation (3) unchanged, but 


results in a change of the following form for the j pertur- 
bation 


L^[<^ 



+ L 
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L[<I>q] = 0 


w 
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Here the operators are understood to be expressed in terms of 
the strained (s,t) coordinates, and the additional operator L2- 

arises specifically from displacement of the strained points. 

In references 3 and 4, specific expressions for L 2 . are provided 

for selected perturbations involving transonic small-disturbance 
and full -potential equation formulations. The essential point, 
however, with regard to perturbation Equation (6) expressed in 
strained coordinates is that it remains valid as before for a 
unit perturbation and independent of . 


In employing the correction method, Equation (6) for the 
j unit perturbation is solved by taking the difference between 
two solutions obtained by the full nonlinear procedure after 
appropriately straining the coordinates. If we designate the 
solutions for some arbitrary dependent flow quantity Q as 
base and calibration , respectively, of the varied 

independent parameter q^ , we have for the predicted flow at 
some new parameter value qj 


where 


M 


Q(x,y) = Q^(s,t) + I e.Q^ (s,t) + 

j=l j 


Qc.Cx^.yp - QoCs.t) 

Qi. = ^ 


N 


X . = s + 
J 


I e.6x.x. (s,t) 

i=l J ^ i 




y = t + j; , (s,t; 

•’ i = l i 


X = s + 


M 0 . 

5=1 j ^ 


- s) 


y = t + 


M e . 

j=l J 


(Continued) 


(7) 
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[Eq. (8) Continued] 


( 8 ) 


In the following section, applications of the correction procedure 
are made to predict surface properties. Also provided is the 
particular form of the straining functions Equation (5] found 
to be most effective in those applications. 


2.3 Application to Surface Properties 

For the current applications, we have employed coordinate 
straining with the correction method to predict distributions of 
surface properties for simultaneous multiple-parameter pertur- 
bations of aerodynamic flows. In that instance where flow 
properties are required along some contour, the strained- 
coordinate solutions can be represented by 

M 

Q(x;c) ~ Qq (s) + I e.Q^ (s) + ... 

j=l J 

M 

X ~ s + e -X-, (s) + . . . 

j=l ^ 

where x is the independent variable measuring distance along 
the contour or a convenient projection of that distance, s is the 
strained coordinate, and Cj a small parameter representing the 

change in one of M flow or geometrical variables which we wish 
to vary simultaneously. 



we require one base and M calibration solutions in which the 
calibration solutions are determined by varying each of the 
M arbitrary independent parameters q^ by some nominal amount 

from the base flow value while keeping the others fixed at 
their base values. 

In this way, the first-order corrections Q]^.(sj can be 
determined as ^ 


(9) 

( 10 ) 
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(113 


(s3 

j 


L_I 

e . 

3 


where . is the calibration solution corresponding to changing 
t h 1 — 

the j parameter to a new value q , x- is the strained coordinate 

^3 

pertaining to the . calibration solution, and e- = q " <ln 

3 ^ ^ 3 ° 3 

represents the change in the q^ parameter from its base flow 

value. If we now desire to keep invariant during the pertur- 
bation process a total of N points corresponding to discontuit ies 
or high-gradient maxima, we can represent the solution by 


M 


Q(x 


= QqCs) + I £-Q;^ (s) 
3 = 1 3 


( 12 ) 


where (s) is given above and 
3 


X - = s 
3 


N 


+ Z i.5x X (s) 
1 = 1 -* 1 


(13) 


X = s + 


N 

^ I e. 6x.x. ^ (s) 
i = l J ■; 


(14) 


£ . = q 
3 ^ • 

J J 


(15) 


e . = q . 
3 ^3 


(16) 


e . 6x . = (x . 

J 1 ^1 


Xi) 


(17) 


£.5x- = -=-^ (x^ - x9) 

3 1 e • 1 1 - 

J £3 J 


(18) 


He 


re £j 5x^ given in Equation (17) represents the displacement 


of the i^^ invariant point in the j calibration solution from 
its base flow location due to the selected change e^ in the q^ 
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parameter given by Equation (15), e.6x. given in Equation (18) 

J ^ - 4 -}, 

represents the predicted displacement of the i invariant point 
from its base flow location due to the desired change in the 

q. parameter given by Equation (16), and x. (s) is a unit-order 
J -*-i 

straining function having the property that 




k = i 

k 5^ i 


(19) 


which assures alignment of the i^^ invariant point between the 
base and calibration solutions. 

In addition to the single condition Equation (19) on the 
straining function, it may be convenient or necessary to impose 
additional conditions at other locations along the contour. For 
example, it is usually necessary to hold invariant the end points 
along the contour, as well as to require that the straining 
vanish in a particular fashion in those locations. All of these 
conditions, however, do not serve to determine the straining 
uniquely. The nonuniqueness of the straining, nevertheless, can 
often be turned to advantage, either by selecting particularly 
simple classes of straining functions or by requiring the strain- 
ing to satisfy further constraints convenient for a particular 
application. 

The fact of nonuniqueness of straining function, however, 
raises a further question of the dependence of the final 
perturbation-predicted result on choice of straining function. 

An initial example of the effect of employing two different 
straining functions for a strongly supercritical flow was 
provided in reference 3, and in reference 1 a detailed examina- 
tion was made of the dependence of perturbation results on 
several classes of different straining functions. Although 
it can be demonstrated (ref. 8) that the final perturbation- 
predicted result obtained when employing strained coordinates 
is formally independent of the particular straining function 
used- -provided that the straining function moves the invariant 
points to the proper locations -- the results of reference 1 
demonstrate that, under certain conditions, particular classes 
of straining functions can induce spurious perturbation results. 
The underlying reason is that, while the perturbation-predicted 
results at and in the vicinity of invariant points are independent 
of the choice of straining function (provided invariant point 
locations are preserved) , some classes of straining functions 
have the undesirable property of producing unwanted straining 
in certain regions removed from the invariant points. The 
correction for this deficiency, which was found in reference 1 
and has proven effective in all case studies undertaken, is to 
employ linear piecewise-continuous straining functions. This 
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both preserves the accuracy of the perturbation results in the 
vicinity of the invariant points, and introduces no excessive 
straining in regions removed from those locations. 

For linear piecewise-continuous straining functions, the 
functional forms of the straining can be compactly written. 


= s + 


i + 1 


o o 

^i+1 ■ ^i 


r C 0 . 

(X^ - xp 


(20) 


o 

^ - ^i 
o o 

- '‘i 


r C O'* 

(Xi^l - x.^i) 




- s) 


H(s - xp 


where H denotes the Heaviside step function. As discussed 
above, it is usually necessary to hold invariant both of the 
end points along the contour in addition to the points corre- 
sponding to discontinuities or high- gradient maxima. Consequently, 
for the results reported here, the array of invariant points 
in the base and calibration solutions have been taken as 


X? = { 0 , x° , X 2 , . . . , X° , 1 } 

X? = {0, , x^ , x^ , 1} 

j j 3 3 


( 21 ) 


where the contour length 
where n is the number of 
contour exclusive of the 


has been normalized to 
invariant points along 
end points . 


unity and 
the blade 
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3. RESULTS 


Because the ultimate utility o£ the perturbation methods 
being developed under this investigation is in engineering 
design or parametric analysis, the two primary objectives of 
the current study were to develop the simultaneous multiple- 
parameter capability of the method and then to examine the 
accuracy and range of validity of the multiple -parameter method 
in situations characteristic of that environment. Toward that 
end, we have tested the method in a series of problems involving 
simultaneous multiple-parameter variations of both flow and 
geometric quantities. As with the testing of the single- 
parameter method (tef . 1) , emphasis was placed on transonic flows 
past both isolated blades and compressor cascades that are 
strongly supercritical and exhibit large surface shock movement 
over the parametric range studied. Additionally, we have 
coupled the mult iple -parameter method with an optimization 
procedure to test the method’s ability to perform in an actual 
design environment. These preliminary case studies of the 
combined perturbation/optimization method actually resulted 
in the most demanding tests of the perturbation method under- 
taken to date for observing its ability to work accurately under 
extreme interpolation/extrapolation conditions. 


3.1 Simultaneous Multiple -Parameter Perturbations 

In Figure 1, we present a comparison of results for the 
simultaneous perturbation of thickness ratio and oncoming Mach 
number of highly- supercritical flows past a series of isolated 
NACA four-digit (OOXX) blade profiles. The base flow chosen 
for these results is at M = 0.820 and t = 0.120, and is 

indicated on both plots shown in Figure 1 as the dashed line. 
Those results were obtained by solving the full -potent ial 
equation based on the finite -difference relaxation approximate- 
factorization method of reference 9. The body-fitted mesh 
employed had 75 points on both upper and lower surfaces. The 
calibration solution selected to account for Mach number changes 
is at M^ = 0.800 and t = 0.120, and is displayed as the dotted 

line in the plot on the left; while the calibration flow 
selected to account for thickness-ratio changes is at M^ = 0.820 

and T = 0.110 and is displayed as the dotted line in the plot 
on the right. The open circles represent the perturbation- 
predicted solution for M = 0.790 and t = 0.115, which is a 

parameter extrapolation in M^ and interpolation in t. Those 

results are meant to be compared with the ’exact’ nonlinear 
results which is indicated as the solid line. We note that the 
indicated results for base, perturbation, and exact nonlinear 
solution in both plots in Figure 1 are the same; the reason 
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for presenting two plots is to indicate clearly the separation 
between the base, the two calibration solutions, and the 
perturbation-predicted result. The straining employed is 
linear piecewise-continuous [see Eq . (20)], with leading and 
trailing edge and shock point held invariant. The shock point 
locations for the base and calibration flows for this example, 
as well as for all the multiple-parameter perturbation results 
presented here, were determined as the point where the pressure 
coefficient passed through critical with comprehensive gradient 

With regard to the results, we note that the comparison 
between the perturbation-predicted and the exact nonlinear 
result is, as in the case of single -parameter perturbation 
of these flows (see Figs. 4 to 6, ref. 1), extremely good, in 
particular in the region of the shock. The multiple-parameter 
perturbation result is able to accurately predict both shock 
location and the critical post-shock expansion behavior. 

Results for the region from the stagnation point to points 
just ahead of the shock are essentially identical to the exact 
nonlinear solution, as are results aft of the post-shock 
region. We note that the particular parameter values of 
(Moo,t) = (0.790 , 0.115) selected for the predi ction solution 
represent reasonably substantial extrapolations and interpo- 
lations from the base and calibration values. Nevertheless, 
the perturbation method is able to treat simultaneous 
parameter variations over this range accurately. 

Figure 2 presents analogous three-parameter perturbation 
results when angle-of -attack variations are included for the 
flows shown in Figure 1. Here, the base flow selected is at 
a = 0.2°, Moo = 0.800 , t = 0.110, and is indicated in all of the 
three plots provided as the dashed line. The calibration flow 
to account for angle-of -attack change is at a = 0.25° at the 
same (Mqo,t) as the base flow, and is displayed as the dotted 
curve in the plot on the upper left. The corresponding 
calibration flow to account for Mach number change is at 
Moo = 0.810 at the same (x,a) as the base flow, and is displayed 
in the upper right plot; while the calibration flow for thick- 
ness-ratio change is at t = 0.115 at the same (Moo, a) as the 
base flow, and shown in the lower plot. The predicted result 
is for parameter values of a = 0.3°, Moo = 0.820, t = 0.100 and 
again represents substantial extrapolations of all three 
parameters, as can be observed in Figure 2 from the relative 
differences between the base and calibration flows. The 
reason for selecting such small angles -of -attack for these 
flows was to preserve the shock wave on the lower surface, and 
thereby create a set of multiple -shock flows which were highly 
sensitive to parameter changes. The comparisons between 
perturbation and exact nonlinear results for this case is again 
extremely good, with the prediction of both the locations of 
the shocks on the upper and lower surface given very well, as 
well as the pressure distributions in the regions immediately 
ahead and behind those shocks. For these results, linear 
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piecewise-continuous straining was employed with the invariant 
points being the lower surface trailing edge, lower surface 
shock, stagnation point, upper surface shock, and upper surface 
trailing edge. 

The final mult iple -parameter perturbation result is 
provided in Figure 3 for a four-parameter perturbation of strongly- 
supercrit ical full -potential flows past a cascade of blades 
having NACA four-digit profiles. The base flow is for an 
oncoming Mach number of Moo = 0.780 , thickness -ratio t = 0.110, 
gap-to-chord ratio t = 3.2, and oncoming inflow angle a = 0.3°, 
and is indicated on the four plots as the dashed line. Those 
results were obtained using the full -potential equation 
finite-area relaxation procedure of reference 10. The four 
calibration solutions to account for changes in the four varied 
parameters are provided in the four plots shown where the 
individual values of the calibration parameter varied are also 
indicated. Thus, the calibration solution for Mach number 
change is at Moo = 0.790 with (T,t,a) at the base flow values, 
and is indicated as the dotted result in the plot at the upper 
left. Corresponding results for the other three parameters 
are shown in the remaining plots. The comparison of the 
predicted and exact nonlinear results are for parameter values 
of Mj„ = 0.785, T = 0.115, t = 3.1, a = 0.4°. This particular 
set of flows was again selected because of the presence of 
multiple - shocks and high sensitivity to parameter change. As 
with the previous results for two- and three simultaneous 
parameter variations, we note that the perturbation predictions 
are once again remarkably accurate. The perturbation method 
is able both to track the location of the upper and lower surface 
shocks, as well as to predict the pressure characteristics in 
the pre-shock and post-shock regions. 


3.2 Preliminary Application of Combined 
Multiple -Parameter Perturbation Method 
With Optimization Procedures For Blade 
Design Applications 

The ultimate utility of the perturbation methods developed 
and evaluated here lies in their application to problems involving 
the high-frequency use of computational codes to determine a 
large number of related nonlinear flow solutions. In order to 
test the capability of the approximation method to work 
effectively in such practical applications, we have combined 
the method with the CONMIN optimization procedure (ref. 11} and 
have then made several preliminary case studies of the combina- 
tion on isolated blade and compressor blade design/optimization 
problems. The objectives of these initial applications were to 
examine the feasibility and potential computational savings of 
the combined approximation/optimization procedure for some 
typical design problems, and to determine the accuracy of the 
perturbation -predicted optimization results. 
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The particular isolated blade design optimization problems 
selected for study involved the alteration of a baseline profile 
shape by adding to the baseline profile a set of shape functions 
according to the relation 

M 

z(x) = (X3 + I (DV. - 1) F^Cxl (22) 

i = l 


where Zq are the ordinates of the baseline profiles, F;j^ are the 
shape functions, and the coefficients DV^ are the design variables 
whose values are determined by the optimization program as a 
result of a search through design-variable solution space to 
achieve a desired design improvement. Here for convenience 
we have chosen the coefficients of to be - 1) rather 

than DVi. The general class of geometric shape functions 
employed here, and which have been found to be successful in 
previous applications involving optimization of supercritical 
airfoil sections (ref. 12) , consists of exponential decay 
functions and sine functions. These are of the general form 
(1 - x) • xP/e^^ and sin (nx^)^, where the exponents p, q, r, 
and n are selected to provide a desired maximum at a particular 
chordwise location. The exponential functions are generally 
employed to provide adjustments near the leading edge, while 
the sine functions are used to provide maximum ordinate changes 
at particular chordwise stations. Illustrations of the chord- 
wise variation of typical members of these shape functions are 
provided in Figure 4, and it can be seen that these functions 
smoothly concentrate ordinate thickness at selected locations. 

For the initial application of the combined perturbation/ 
optimization method, we have examined subcritical flow at 
Moo = 0.10 and a = 5° past a modified NACA 64A007 profile 
involving the nine profile shape functions 


F- = 6(1 - x)xPi/e^li^ i = 1,2 
r- 2 

= sin (ttx 1) i - 3,9 


(23) 


where (p^,q^) = (0.5,15), (p 2 ,q 2 ) = (0.25,10), and r^ = (0.37, 

0.50, 0.66, 0.87, 1.16, 1.61, 2.41). The exponential functions 
achieve their maxima within 5°g of chord, while those for the 
sine functions are at (15%, 25%, 35%, 45%, 55%, 65%, 75%) of 
chord . 

A strategy which has proved convenient for performing 
optimization studies involving aerodynamic performance 
parameters (ref. 12) has been to recontour the profile shape 
so as to tailor the surface pressure distribution to conform to 
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a desired distribution. This type of objective provides local 
control over the basic aerodynamic surface flow property of 
importance, and provides a means of attempting to achieve aft 
pressure gradients sufficiently weak to avoid separation. An 
important corollary advantage of using such an objective is that 
viscous separation can be minimized. This allows use of an 
inviscid aerodynamic flow solver in the optimization process 
rather than a much more computationally- expens ive viscous solver, 
and assures that the optimization result thus obtained at the 
inviscid level is representative of the actual flow. 

Consequently, for this initial case study the overall 
performance objective was, through modification of the surface 
contour, to tailor the pressure distribution along a portion of 
the upper surface so as to conform to a desired distribution. 

In particular, it was desired to minimize both the peaky 
behavior near the leading edge and the compressive gradient 
on the aft portion of the upper surface which existed at Moo = 0.10 
and a = 5° on the NACA 64A007 baseline profile. This is illus- 
trated schematically in Figure 5. The objective function was 
taken as the minimization of the mean squared error between the 
predicted and desired surface pressure distribution, i.e.. 
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where K represents the number of chordwise locations where 
desired and calculated surface pressures are compared. 

Recall that in order to initiate the perturbation proce- 
dure in situations involving the simultaneous variation of M 
individual parameters from a baseline point, a matrix of M 
calibration solutions are required, each representing the 
solution change for a separate variation of each of the M 
parameters from its baseline value. Because optimum, or some- 
times even typical, stepsizes for a particular optimization 
problem would not generally be known a priori, one of the 
primary goals of these initial studies was the demonstration 
that the perturbation method was capable of working effectively 
even under severe conditions imposed by a poorly- selected 
initial calibration solution matrix. This was accomplished by 
examining the sensitivity and accuracy of the perturbation 
predicted optimization results as a function of the initial 
design variable stepsizes of the calibration solution matrix. 

The overall strategy, then, consisted of: (1) employing the 

perturbation method, based on some initial matrix of nonlinear 
aerodynamic solutions determined by an independent variation of 
each design variable, to provide all of the subsequent nonlinear 
aerodynamic solutions required by the optimization method for 
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searches through design space, and (2) comparing the final 
perturbation-predicted optimization results for final design 
variable and objective function values with those determined 
by using the full nonlinear aerodynamic solver throughout. 

Figure 6 shows the results of such a sensitivity study, and 
indicates that even under extreme test case conditions caused by 
deliberately-selected poor choices of design variable stepsizes, 
the perturbation method performs exceptionally well, never 
breaking down or yielding spurious results. Indicated on the 
plots are the final optimized design variable values after 
five search cycles as predicted by the perturbation method (#) 
for four different choices of the initial stepsize for the 
design variables, i.e., <SDVj^ = (0.05 , 0.02 , 0.01 , 0.001). Also 
shown are the corresponding final design variable values 
predicted when employing the nonlinear aerodynamic solver 
(ref. 13) throughout (O). As can be seen, for the extreme 
interpolation case = 0.05, except for design variables 3 

and 5, the perturbation prediction compares very favorably with 
the full nonlinear result. For 6DVj^ = 0.02, a more reasonable 
stepsize choice, the comparison of the perturbation result is 
quite good for all the design variables, while for 6DVj^ = 0.01, 
the perturbation prediction is essentially identical to the full 
nonlinear aerodynamic result. As a final illustration of the 
behavior of the perturbation method under extreme extrapolation 
conditions, the lower right-hand plot displays the perturbation 
predictions for 6DVj^ = 0.001. We note that for several of the 
design variables, the extrapolation range is of the order of 25 
times the initial stepsize; yet, the perturbation predictions 
are quite reasonable and not spurious, which is remarkable and 
indicative of the robustness of the procedure. 

Finally, we point out that all four of the perturbation- 
predicted results illustrated in Figure 6 are satisfactory in 
terms of the final objective function value obtained. These 
values are illustrated on the right of each of the plots. 

Provided for comparison are the initial (O Initial) and final 
(O Final) values obtained when using the nonlinear aerodynamic 
solver throughout. The value of the objective function evaluated 
at the final design variable point when using the perturbation 
method is indicated by the solid circle (•). However, the 
objective function value of real interest is the result indicated 
by the solid square (■) which represents that obtained by 
running the nonlinear aerodynamic solver at the perturbation- 
predicted final design point, and then using that solution to 
evaluate the objective function. This provides the overall 
ultimate check of the perturbation-predicted result. As can 
be seen from Figure 6, those results lie essentially on top of 
the final objective function result (O Final) obtained when 
using the nonlinear aerodynamic solver throughout. 

The computational savings attained for this application 
are shown in Figure 7. Here, a comparison of the computational 
work versus reduction in objective function per optimization 
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cycle is provided when using the perturbation procedure (#) and 
when not using it but employing the nonlinear aerodynamic code 
(O) for each flow solution required by the optimizer. As 
can be seen, the computational time required for both the 
perturbation method and when using the full nonlinear aerody- 
namic solver throughout are the same for the first cycle, 
since both require a matrix of M + 1 (M=9 for this example) 
nonlinear aerodynamic solutions. After that, the perturbation- 
predicted results required essentially no computational time 
for cycles 2 through 5, and then a slight amount for the one 
additional call to the aerodynamic solver for the final check 
calculation (■)• The reduction in the ratio of final to 
initial objective function is OBJ/ (OBJ) = 0.22 and required 
approximately 20 CPU seconds on the CDC 7600. In comparison, 
the result when not employing the perturbation method required 
approximately 80 CPU seconds for the same reduction in 
objective function^ indicating the perturbation method is able 
to save 75% of the computational work in this example. 

Similar testing of the performance of the perturbation 
method for supercritical situations has also been carried out 
and has demonstrated a corresponding capability and potential 
computational savings. Because of the greater sensitivity of 
these shocked flows, two separate case studies were carried out 
in depth. Both of these involved recontouring of the upper 
surface of a NACA 0015 profile operating at the supercritical 
conditions of Moo = 0.55 and a = 6.7°, and employed four design 
variables related to the shape functions . 

q . 

= sin(7Tx i = 1,4 (25) 

with q. = (0.301, 0.431, 0.576, 0.756). These functions have 
maxima^at (10%, 20%, 30%, 40%) of chord. For the first of 
these supercritical studies, the objective function was chosen 
to be the drag coefficient squared, i.e. 

OBJ = CD^ 

For this problem, an evaluation of the accuracy of the 
perturbation-predicted optimization results as a function of 
initial calibration solution matrix was also made. The results 
of this study are provided in Figure 8, which displays the 
results of the perturbation-predicted final design variables 
(#) after eight search cycles for three different choices of 
the initial stepsize for the design variables, i.e., = 

(0.001, 0.002, 0.004). Also shown are the corresponding final 
design variable values predicted when not using the perturbation 
method but . employing the nonlinear aerodynamic solver (ref. 13) 
throughout ( O ) . 


17 



Because drag minimization, particularly at supercritical 
conditions, is an extremely sensitive optimization problem, 
this study provides one of the most severe tests of the 
perturbation procedure in a design optimization environment. 

This is so because the accurate prediction of drag at super- 
critical speeds depends almost entirely on the resolution of 
the flow behavior at and in the vicinity of the shock waves 
present on the surface of the profile. Hence, what is ultimately 
under evaluation in this example case study is the ability of 
the perturbation method to predict, under extreme extrapolation 
conditions and with simultaneous multiple-parameter perturbations, 
the location and strength of all surface shock waves and the 
flow behavior in the pre- and post-shock regions. 

The most important results to emerge from the calculations 
involved in this case study were the discovery of a particular 
deficiency of the perturbation method in this regard, and the 
subsequent development of the means to improve the accuracy of 
the perturbation predictions in shock regions and other high 
gradient regions under extreme extrapolation conditions. The 
improvement in the basic procedure developed to meet these 
requirements consists of employing additional invariant points 
in those high gradient locations. For example, it was found 
that by characterizing a shock which has a post-shock expansion 
region, as sketched below. 



with five invariant points - -which correspond tcK pre-shock 

minimum pressure, © maximum gradient point, (3) post-shock 
maximum pressu:^, (3) post-shock minimum expansion pre^ure 
l^ation, and © point of inflection between points (?) and 
( 4 ) --rather than just the one point corresponding to the 
critical pressure location, which was standardly done in the 
past; that significantly improved perturbation results are 
obtained in the shock region for extreme solution extrapolations. 
This five invariant point characterization of the shock has been 
employed in determining the results of the two supercritical 
case studies reported here. 
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With regard to the drag minimization results indicated in 
Figure 8, we note that when selecting an initial design variable 
stepsize o£ 6DVj^ = 0.001, the agreement between the final 
perturbation-predicted design variable values and the exact 
nonlinear result is reasonable. The solution extrapolation 
indicated for design variables 2 and 3 are of the order of 5 
times the initial design variable stepsize, but the perturbation 
method does not break down or provide spurious results for these 
strongly supercritical flows. We note that the optimization 
procedure working with the perturbation method is able to drive 
the perturbation-predicted drag P) to essentially zero. 
Although the final check of the perturbation-predicted design 
using the nonlinear solver (■ A) indicates a drag coefficient of 
0.005, that represents nevertheless almost a factor of 3 in 
drag reduction from the baseline configuration. 

When a somewhat more reasonable initial design variable 
stepsize of (SDVjL = 0.002 is used, which reduces the solution 
extrapolations over that for 6DVj^ = 0.001, we note that the 
perturbation results closely approach those obtained when 
using the nonlinear aerodynamic solver throughout. The final 
drag value based on the perturbation design (BA) is also 
improved, displaying over a factor of 5 in drag reduction. The 
final perturbation result obtained for an initial stepsize 
6DVi = 0.004 and illustrated in the bottom plot displays almost 
exact agreement with the full nonlinear aerodynamic result. 

This particular choice of stepsize is the most reasonable of 
the three shown, since the design variable solution range to 
be searched by the optimizer requires moderate interpolation/ 
extrapolations from the calibration solution matrix. The final 
drag value (BA) indicates over a factor of 5 drag reduction, 
and is very close to the final full nonlinear aerodynamic 
(G Final] result. 

The computational savings obtained for this application 
are indicated in Figure 9. Here a comparison of the computa- 
tional work versus reduction in objective function per 
optimization search cycle is given for eight search cycles 
when using the approximation procedure (•) and for four cycles 
when using instead the full aerodynamic solver (O). We note 
for this example that when using the perturbation method, the 
optimization procedure requires approximately 40 CPU seconds 
for search cycles 2 through 8. This is in contrast to the 
essentially zero time needed for search cycles in the pressure 
distribution tailoring case study presented in Figure 7. The 
reason for the additional time in this minimization case study 
is that in order to evaluate the drag at the points in the 
design solution space required by the optimizer, additional 
geometry calculations are needed to determine the new surface 
slopes. Consequently, even though the perturbation procedure 
provides the means to determine the new pressure distributions 
at essentially no computational cost, because surface pressure 
integrations are required to evaluate the drag, calls to the 
geometry portion of the aerodynamic code are necessary to 


determine the new profile shapes. These geometry calculations 
cause the CPU time increase indicated. Finally, another 20 
CPU seconds are needed for the final aerodynamic check solution 
(■), totaling 120 CPU seconds. That total represents just 
under a 60% computational savings of the 4 search cycle CPU time 
required when not employing the perturbation method. 

The final supercritical case study presented in Figure 10 
involved tailoring of the surface pressure distribution as done 
in the previous subcritical case, but this time specifically 
concentrating on alleviating the sharp gradient at the shock. 

It was specified in the following fashion. Using the pressure 
distribution from the optimized profile of the previous drag 
minimization case study, we returned to the original NACA 0015 
baseline profile. Then, with that previous pressure distribution 
as the desired objective, by using the optimization method we 
attempted to reproduce the same optimized profile that resulted 
from the drag minimization problem. What resulted from these 
calculations was a demonstration that the perturbation method 
is capable in certain cases, of not only providing an enormous 
savings in computational cost but also an improved optimization 
result. In Figure 10, we provide the results of the sensitivity 
study of the perturbation method as a function of the initial 
design variable stepsize of the calibration solution matrix. 

In the plot on the left, the comparison is provided of the 
perturbation-predicted final design variables (•) based on 
an initial design variable stepsize 6DVj^ = 0.001 together with 
that obtained when using the full nonlinear aerodynamic solver 
throughout C^). Also shown is the ’optimum’ full aerodynamic 
result (O) obtained from the drag minimization case study which 
is the result that is sought to be reproduced. What the results 
indicate is that the prediction obtained by not using the 
perturbation method is inferior to that obtained when employing 
it. Also, we note the large extrapolations involved in this 
perturbation result, i.e., about 4 times or 400% of the initial 
stepsize for design variables 1 and 2. These results are 
emphasized even more in the plot on the right which shows the 
corresponding perturbation results when using a more reasonable 
initial design variable stepsize of 6DVj^ = 0.002 to define the 
calibration solution matrix. In this instance, the perturbation 
prediction is essentially identical to the optimum result. We 
note that comparisons of the objective function reduction, shown 
on the scales to the right of the plots, display a reduction in 
objective function to almost zero for the perturbation result, 
while only an essentially 50% reduction for the full aerodynamic 
result . 

This reduction in objective function is emphasized in 
Figure 11, which displays the comparison of computational work 
and objective function reduction per optimization cycle for 
the perturbation procedure and the full nonlinear aerodynamic 
result for this case study. Here, we see clearly that the 
perturbation method is able to drive the objective function to 


20 



essentially zero, while the full aerodynamic procedure becomes 
fixed in a local minimum and is only able to reduce the 
objective function to 50% of its initial value. If we had 
carried the full aerodynamic result to eight optimization cycles, 
as was done with the perturbation result, the time savings would 
have been over an order of magnitude. 

This result demonstrates that it is possible in certain 
instances for the perturbation method to provide not only an 
enormous savings in computational cost but also an improved 
optimization result. The latter is undoubtedly accomplished 
by the selection of a reasonable initial calibration solution 
matrix which permits the optimization method an enhanced 
rather than local view of the design solution space, thereby 
avoiding shallow local minimas in favor of a more global minima. 

The final application, which has been carried out as far 
as possible in this phase, was directed toward laying the 
foundations of a practical turbomachinery blade design/ 
optimization procedure coupled with the perturbation method. 

This preliminary version is based on TSONIC blade- to-blade 
solutions (ref. 14) with generalized circular-arc blade 
geometry routines BLADE (ref. 15) describing the blade profile, 
and employs the COPES/CONMIN optimization procedure (ref. 16). 

The initial case study for the combined procedure involved, 
as a design objective, the minimization of the velocity diffusion 
on the blade suction surface, q^ax, suction’ 


OBJ 


*^max, suction 
‘^avg , exit 


(27) 


where q is the average exit velocity in the freestream. 

^avg , exit 

Six design variables were employed and correspond to the 
following blade geometry parameters used to characterize NASA 
circular arc blade section profiles (ref. 15): blade outlet 

camber angle, transition location between fore and aft circular 
arc sections, maximum thickness location, inlet to outlet 
turning rate ratio, maximum thickness, and radius of the lead- 
ing edge circle. During the optimization process, each of the 
design variables was constrained to remain within certain 
prescribed bounds in order to prevent a physically-unrealistic 
blade design from occurring. Furthermore, several active side 
constraints were additionally imposed both to insure design of 
a physically realistic blade and also to achieve certain 
desirable flow characteristics on the blade. The active side 
constraints employed were: (1) maintenance of nonzero local 

blade thickness, (2) maintenance of low velocity diffusion on 
the blade pressure surface, and (3) trailing edge closure. 

These constraints were enforced by employing bounding criteria 
on the following quantities according to: 
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thick (x) 


< 10 . 


0 . < 



(28) 


0. < press ^ 

^min, press 


(29) 


1 0 < 2 , suction ^ITE-2, press 

25 . 


< 1.0 


(30) 


Here, thick (x) is the local blade thickness, is the diameter 

of the trailing edge circle, is the maximum blade 

IU3.X j p T© S S 

pressure surface velocity over the front half of the blade, 

is the minimum blade pressure surface velocity 
min, press ^ ^ 

over the last two- thirds of the blade, and (qr^rT- o 

> ''^ITE- 2 , suction, 

^ITE-2 press^ third last surface velocities on the grid 

near the trailing edge on the suction and pressure surface, 
respectively . 


We have successfully completed a preliminary series of 
calculations of the new combined PERTURB/TSONIC/BLADE/COPES/ 
CONMIN procedure in which the accuracy and sensitivity of the 
perturbation method was tested as a function of choice of the 
initial calibration solution matrix. The initial or base values 


of the design variables for the baseline blade 
upper and lower bounds of the design variables 
specified for this test problem were: 

profile , 
that were 

and the 

Design 

Variable 

Number 

Description 

Lower 

Bound 

Upper 

Bound 

Initial 

Value 

1 

Outlet blade camber angle-KOCR 

-15.0° 

0.0° 

-10.0° 

2 

Transition loc./chord-T 

0.20 

0.60 

0.25 

3 

Max. thick, loc . /chord- ZM 

0.20 

0.55 

0.45 

4 

Inlet/outlet turning/chord-P 

0.50 

4.00 

1.50 

5 

Max, thick . /chord-TMX 

0.03 

0.10 

0.05 

6 

L.E. rad . /chord-THLE 

0.003 

0.012 

0.005 

The results of these calculations are summarized in 

Table 1. 


There we have provided comparisons of the final design variables 
and objective function predicted when employing full nonlinear 
TSONIC solutions throughout the optimization process with 
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corresponding results when using the perturbation method. For 
the perturbation results, six different choices of the cali- 
bration solution matrix were made and are noted in the table. 

All the results represent converged solutions, with each 
calculation employing 10 optimization search cycles or less 
if no change in objective function should occur in three 
successive iterations. 

Two separate results are presented for the case when the 
perturbation method is not employed and TSONIC solutions are 
used throughout the optimization search process. These were 
obtained on the Ames Research Center CDC-7600 and Lewis Research 
Center IBM 3033 employing the same TSONIC/BLADE/COPES/CONMIN 
code. The differences in the two final design results provides 
an indication of the sensitivity of this particular optimization 
problem with regard to choice of objective function, since the 
difference between these two results is due solely to the 
number of significant figures maintained in the respective 
calculations, i.e., 8 for the IBM 3033 and 14 for the CDC 7600. 

As is evident, there clearly is a sensitivity to search direction 
and final design result in this problem. We note that of the 
six design variables, the two which agree between the two 
results (ZMjTMX) have reached a limit boundary. The others 
have all trended in the same direction from the baseline value, 
but have reached somewhat different values at the final design 
point. We note that the objective functions have reached 
essentially the same level, indicating, as is often the case 
in such optimization problems, the existence of multiple local 
minimums . 

The analogous results obtained employing the perturbation 
method with various choices of calibration solution matrix 
are exhibited as cases 1 to 6. We note that, with the exception 
of the design variable ZM which consistently moves to its 
upper bound regardless of the choice of calibration solution 
matrix, the final design variable predictions via the perturbation 
method basically exhibit a behavior quite similar to that 
obtained from the full nonlinear result when TSONIC solutions 
are used throughout. That is, the optimizer drives the design 
variables in generally the same direction from the baseline 
values as the full nonlinear result, but to somewhat different 
values. The exception to this is noted in cases 4 and 5, 
where the calibration solution matrices were selected such that 
TSONIC solution extrapolations rather than interpolations were 
involved for all or most of the design variables during the 
search procedure. In those two cases, we note further that 
the final objective function values obtained are somewhat 
inferior to those obtained in Cases 1, 2, 3, and 6, where 
solution interpolations were primarily involved. For those four 
cases, the final objective function result is almost identical 
to that obtained by the full nonlinear result. 
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The computational time needed to obtain the perturbation 
results in Cases 1 to 6 were 76-78 secs, of CDC 7600 CPU time 
per case. The benchmark full nonlinear CDC 7600 result shown 
in Table 1 required 644 secs. Thus, the perturbation method 
provides a savings of (644-78)/644 = 88% of the computational 
time for this example. 

The primary conclusions to be drawn for this preliminary 
study are that the perturbation method can work effectively 
even for sensitively-defined optimization problems such as 
this and provide both meaningful final design results and 
large computational savings over not using the method. The 
choice of objective function such as was made for this case 
study, i.e., a point quantity located in a high- gradient 
region, requires that a reasonably good choice be made of the 
calibration solution matrix. This is so since if large solution 
extrapolations are required by the optimizer for minimization 
searches through design space, the final design result will 
most likely be less improved than would otherwise result if only 
modest solution interpolations/extrapolations were involved. 
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4. CONCLUSIONS AND RECOMMENDATIONS 


An investigation was conducted to continue the development 
o£ perturbation procedures and associated computational codes 
for rapidly determining approximations to nonlinear flow 
solutions. The ultimate purpose is to establish a method for 
minimizing computational requirements associated with para- 
metric design studies of transonic flows in turbomachines. 

The procedures being developed employ unit perturbations, 
determined from two or more nonlinear ’base' solutions which 
differ from one another by a nominal change in some geometry 
or flow parameter, to predict a family of related nonlinear 
solutions which can be either continuous or discontinuous. The 
results reported here relate to the extension of the previously- 
developed successful method for s ingle -parameter perturbations: 

(1) to simultaneous multiple-parameter perturbations, and (2) 
to the preliminary application of the multiple-parameter procedure 
in combination with an optimization procedure and applied to 
blade design problems. 

Calculations based on ful 1 -potential nonlinear solutions 
have been carried out to establish the accuracy and range of 
validity of the mult iple -parameter capability. These involved 
flows past both isolated blades and compressor cascades involv- 
ing simultaneous changes in both flow and geometric parameters, 
with attention focused on strongly supercritical situations 
involving large surface shock movements over the parameter 
ranges studied. Preliminary applications of the multiple- 
parameter perturbation method coupled with an optimization 
procedure were made for blade design problems in order to 
examine the capability of the method to produce accurate results 
in a typical design environment, and also to evaluate its 
potential for computational savings. Both subcritical and 
supercritical case studies were carried out involving multiple- 
design variables. Sensitivity studies were also performed 
to examine the accuracy dependence of the perturbation method 
on the choice of the initial calibration solution matrix. 

Comparisons of the multiple-parameter perturbation results 
with the corresponding 'exact' nonlinear solutions display 
remarkable accuracy and indicate a large range of validity, in 
direct correspondence with previous results for single-parameter 
perturbations. The preliminary case studies of the multiple- 
parameter perturbation method combined with optimization 
procedures have clearly demonstrated the capability of the method 
to work accurately in a design environment where large solution 
extrapolations are often necessary, and have also established 
the methods' potential for reducing the computational work 
required in such applications by nearly an order of magnitude. 

The sensitivity studies indicate that for certain subcritical 
applications, the perturbation method is able to work quite 
accurately and effectively in spite of poor choices of the 
initial calibration solution matrix which require large extra- 
polations or interpolations. For supercritical flows, the 
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initial calibration matrix choice is more important, but when 
employing multiple invariant point clustering about high- 
gradient locations (shocks, stagnation points, etc.), the 
perturbation method predictions can nevertheless maintain 
high accuracy in certain supercritical situations for extra- 
polations as large as 4 to 5 times the parameter separation 
between base and calibration solutions. 

Based on these results, we conclude that perturbation 
methods formulated on these ideas are both accurate and clearly 
workable in design environments, and can provide the means 
for substantially reducing the computational work required in 
such applications. We suggest the development of the combined 
multiple-parameter perturbation procedure and optimization 
methods into a robust design tool for subcritical and super- 
critical turbomachinery blade design. 
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APPENDIX A 

USER'S MANUAL FOR COMPUTER PROGRAM PERTURB 


A . 1 INTRODUCTION 

The purpose of this appendix is to describe the operation 
of the computer code which was developed in conjunction with 
the theoretical work presented in this report, and to provide 
sufficient detail to permit convenient use and change of the 
program. The program computes and plots an arbitrary flow vari- 
able on a contour surface by employing the strained-coordinate 
perturbation method discussed in the main text. 

A description of the general operating procedure of the 
program is given, together with complete description of both 
input and output. The program is written in FORTRAN IV and has 
been developed on the Ames Research Center CDC 7600 computer 
facility. Consequently, the plot package included in this 
version refers to system routines at that facility. In general, 
the plotting software must be supplied by the user according to 
the requirements of his operating system. This can be accom- 
plished directly by replacing or modifying the subroutine DRVPLT . 
Typical program run times for flows involving approximately 150 
surface mesh locations are 1 to 3 CPU seconds. The storage 
requirements are lOSKg for small core memory and no large 
core memory. 


A. 2 PROGRAM DESCRIPTION 

The program calculates both continuous and discontinuous 
nonlinear perturbaton solutions which represent a multiple- 
parameter change in either geometry or flow conditions by 
employing a strained-coordinate procedure. The method utilizes 
unit perturbations, determined for each parameter from a previ- 
ously-calculated common base solution and a calibration solu- 
tion displaced from it by some reasonable change in the relevant 
parameter, to predict new nonlinear solutions over a range of 
parameter variation. 

This current |version of the procedure is configured to 
predict and plot an arbitrary flow variable (e.g., pressure coef- 
ficient) on the surface of a blade or airfoil, and can account 
for the motion of : 

1. one or more critical points (shock points) , 
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2. a stagnation point, 

3. a maximum-suction-pressure point, 

or simultaneously for any combination of these. 

The program is also configured to compare the perturbation- 
predicted solutions with the corresponding ’exact’ solutions 
obtained by employing the same ’expensive' computational proce- 
dure used to determine the base and calibration solutions. 

The coordinate straining employed is piecewise linear 
with the end points and up to six interior points held invariant. 
At the option of the user, these additional interior points may 
be arbitrarily preselected, or chosen from among the minimum, 
maximum, and critical points automatically located by the 
program itself. 

Critical or shock points are located on the basis of a user- 
supplied statement function defining the critical value of the 
dependent variable as a function of some single flow variable. 

The program default is with dependent variable y defined as 
pressure coefficient, with the independent variable being Mach 
number. In this case, the critical value is defined as 


r = c* = 

crit p 
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where y is the ratio of specific heats. If instead of surface 

pressure coefficient, the surface velocity distribution were 

used, then the value of y would be given by 
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Data for base, calibration, and comparison solutions (if 
available) are input as an array x(I) of coordinates a corre- 
sponding array y(I) giving the dependent variable at each coordi- 
nate location, where 1 < I < N and N < 200. 
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I 


Oncoming 

flow 



The leading edge is at x = 0; the data are read in beginning on 
the lower surface at the point farthest from the leading edge 
and proceeding clockwise around the surface as shown in the 
sketch. Data for the different solutions need not correspond to 
identical locations on the surface, except for the initial and 
final points, i.e., x(l) and x(N) must be the same for all cases. 
The program normalizes the x coordinates (0 ^ x ;< 1) such that 
X = 0 corresponds to I = 1 and x = 1 to I = N. 

The base and calibration solutions are searched for minimum, 
maximum, and critical points, e.g.. 
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#2 


(normalized ) 


Direction of 
flow 


Note that the sign of dy/dx in physical coordinates is used in 
determining the critical points. For example, both critical 
points indicated on the above figure correspond to dy/dx < 0 in 
physical coordinates, since at point #1 the physical coordinate 
increases in the direction from right to left, whereas at point 
#2 it increases from left to right. 

The points to be held invariant in straining are either 
selected from among those (1) located by the program or (2) 
individually specified by the user, after which the unit coordi- 
nate straining and unit perturbation are computed. 

Data for the test cases is then read in and nonlinear 
perturbation solutions constructed from the unit perturbation. 
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A 3. PROGRAM FLOW CHART 
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A. 4 DICTIONARY OF INPUT VARIABLES 


A 


B 


LCHEK 


LECHO 


LOCO (I) 


LOCI (I) 


LPLOT 


LSELCT(I) 


Scaling parameter in straining procedure. A = -x(l), 
where x(l) is location of first data point on lower 
surface (see PROGRAT4 DESCRIPTION) . 

Scaling parameter in straining procedure. B = x(N ) , 
where x(n) is location of last data point on upper 
surface (see PROGRAM DESCRIPTION) . 

Specifies whether or not perturbation solution is to 
be checked against an exact comparison solution. A 
printer plot is made in either case. 

LCHEK = 0 ... no comparison 
LCHEK =1 ... comparison 

Controls whether or not input deck is printed. 

LECHO = 0 ... no print 
LECHO =1 ... print 

Array of length 6 of which NSELCT elements are read 
in. Specifies subscripts of those user-specified 
base flow points which are to be held invariant; 
operational only when LSPEC = 1. 

Array of length 6 of which NSELCT elements are read 
in. Specifies subscripts of those user-specified 
points in the Kth calibration solution which are to 
be held invariant; operational only when LSPEC = 1. 

Specifies whether or not an additional plot by a 
peripheral device is to be made. Software must be 
supplied by user in subroutine DRVPLT . 

LPLOT =0 ... No peripheral plot 
LPLOT =1 ... Peripheral plot 

Array of length 6 of which NSELCT elements are read 
in; operational only when LSPEC = 0, and specifies 
nature of points to be held invariant according to 
the code: 

1 ... minimum point held invariant 

2 ... maximum point held invariant 

3 ... 1st critical point held invariant 

4 ... 2nd critical point held invariant 

5 ... 3rd critical point held invariant 

6 ... 4th critical point held invariant 

Note that critical point ordering is determined 
from order of occurence starting at the lower 
surface at the point furthest from the leading 
edge and proceeding clock-wise around the surface 
(see PROGRAM DESCRIPTION) . 
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LSPEC 

LUNIT 

MO, Ml, M2 
N 

NCASE 

NPARAM 

NSELCT 

PARNAM (K) 
QO(K) 

Q1 

Q2 (K) 
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Note that the code numbers can be assigned in any 
order , e . g , , 

LSELCT(l) = 1 LSELCT(l) = 4 

LSELCT(2) = 3 and LSELCT(2).= 1 

LSELCT(3) = 4 LSELCT(3) = 3 

are equivalent, both corresponding to NSELCT - 3, 
with the minimum, and first and second critical 
points held invariant. 

Controls how invariant points in straining are 
specified. 

LSPEC = 0 . . . invariant points selected from 
among those located by the 
program using the array LSELCT(l) 

LSPEC =1 ... invariant points preselected by 
user, using the arrays LOCO (I), 
LOCI (I) 

Controls whether or not unit coordinate strainings 
and unit perturbation ( s) are printed. 

LUNIT = 0 ... no output 
LUNIT =1 ... output 

Oncoming Mach numbers in base, calibration and 
computed perturbation solutions. 

Number of locations for which data are input for 
base, calibration, and comparison solutions. 

Number of cases for which perturbation solutions are 
to be computed . 

Number of parameters perturbed. 

Number of points (in addition to end points) to be 
held invariant in straining; note: 1 £ NSELCT £ 6. 

Array of 8~character strings which identify the 
parameters varied. NPARAM elements of the array 
are read in . 

Array of length 8 giving values of perturbation 
parameters in base solution. NPARAI'4 elements of the 
array are read in. 

Value of Kth perturbation parameter in Kth 
calibration solution. 

Array of length 8 giving values of the perturbation 
parameters in the solution to be computed. NPARA^^ 
elements of the array are read . in. 



TITLE 


Character string of length 80; identifies job and 
is printed as headline on first page of output. 
First nine characters are printed in upper-right- 
corner of banner page, and in upper-left corner 
of summary page. 

VNAM Character string of length 2 which symbolizes 

dependent variable, e.g., "CP" for pressure coeffi- 
cient . 

XBASE (I) ,XCALB (I) ,XCHEK (I) . . . 

Arrays of surface coordinates in base, calibration, 
and comparison solutions. 

YBASE(I) ,YCALB(I) ,YCHEK(I) 

Arrays of dependent variables in base, calibration, 
and comparison solutions. 


A. 5 PREPARATION OF INPUT DATA 
A. 5.1 Description of Input 

Item 1 One card, identifies job-printed as headline on first 
page of output. First nine characters are printed in 
upper-right corner of banner page, and in upper-left 
corner of summary page. 

Item 2 One card, containing the parameter, LECHO 

Item 3 One card, containing the parameters N, NCASE, NPARAM, 

NSELCT, LSPEC, LUNIT , LCHEK , LPLOT . 

Item 4 One card, containing the parameter array (LSELCT(I) , 
1=1, NSELCT) . This item omitted if LSPEC = 1 

Item 5 One card, containing the character string VNAM. 

Item 6 One. card, containing the character strings, 

PARNAM(K) , K = 1, NPARAM 

Item 7 One card, containing the scaling parameters A and B. 

Item 8 One card, containing the parameter MO 

Item 9 One card, containing the parameter array Q0(K), k = 1 

NPARAM ' 
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Item 10 One set of J cards, where J = 1 + INT(N/8) , containing 
data for the base flow coordinates XBASE (I) , I = 1, N 

Item 11 One set of J cards, J as above, containing data for the 
dependent variable in the base solution, YBASE(I), 

I = 1, N 

Item 12 One card, containing the parameter array LOCO (I) ,1=1, 
NSELCT. This item omitted if LSPEC = 0 

Item 13 One card, containing the parameters Ml, Q1 

Item 14 One set of J cards , J as above , containing data for the 
coordinates in the Kth calibration solution, XCAL(I), 

I = 1, N 

Item 15 One set of J cards, J as above, containing data for the 
dependent variable in the Kth calibration solution, 

YCALB (I) , I + 1, N 

Item 16 One card, containing the parameter array LOC(I), 1=1, 
NSELCT. This item omitted if LSPEC = 0 

Item 17 One card, containing the parameter MC 

Item 18 One card, containing the parameter array Q2 (k) , K = 1 

NPARAM 

Item 19 One set of J cards, J as above, containing data for 

the coordinates in the comparison solution, XCHEK(I), 

I = 1, N. This item is required only when LCHEK = 1. 

Item 20 One set of J cards, J as above, containing data for 

the dependent variables in the comparison solution, 
YCHEK(I) , I = 1, N. This item is required only when 
LCHEK = 1. 

NOTE: One set of items 13 through 16 is required for each of 

the NPARAM calibration solutions. 

One set of items 17 through 20 is required for each of 
the NCASE solutions to be computed . 
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A. 5. 2 Format of Input Data 

Item no. 1: 

1 card (8A10) 

Variable 


Title 

Card column 

in 

2nl Tnl 4ol sol . fin 1 8ft 

Format Type 


a 


Item no. 2 ; 

1 card (15) 

Variable 

LECHO 

7 

Card column 

5 

\ 

Format type 

I 

i 


Item no. 3 : 
Variable 
Card Column 

Format Type 


Item no. 4 ; 
Variable 
Card column 

Format Typr 


1 card (1615) 


N 

NCASE 

NPARAM 

NSELCT 

LSPEC 

LUNIT 




5 

10 

15 

20 

25 

30 

35 

40 

VI 

I 

I 

I 

I 

I 

I 

I 

I 

■ 

(LSPEC-Ql: 1 card .f 16151 .^T..qRT.rT fTJSETr-Tl 


LSELCT(l) 





HHHHI 

miHHii 


5 

10 

15 

20 

25 

30 

35 

n 

I 

I 

I 

I 

I 

I 

I 

n 


Item no. 5 : 1 card (2Al) 


Variable 

VNAM 

7 

Card col\imn 

2 


Format type 

A 

i 


Item no. 6 
Variable 

Card column 
Format type 


PARNAM ( 1 ) 

PARNAM (2 



—F 



ft 

16 

24 

32 

40 

48 

56 

A 

A 

A 


HWi 

HBHH 

A 










































Item no. 7. 1 card (8F10.6) 

Variable A b 

Card column ^ 

Format type F F 


> 

CD 

7 

10 20 


F F 

1 


Item no. 8. 1 card (8F10.6) 

Varible i wu > 

Card column 10 \ 

Format type F I 


Item no. 9. 1 card (8F10.6) 

Variable QQ ( 1 

Card colvimn 101 ^ 

Format type F | F 


QO(NPARAM) 


40 1 5 


F F 


Item no. 10. J cards, J»1 +INT (n/ 8) , 8 values per card (8F10.6) 
Variable [xBASE (1)1 XBASE (2)1 XBASE (3 1 

Card column 10 2^ iQ 4^ ^ 

Format type F F F F F 

Item no. 11. J cards, J as above, 8 values per card (0F1O.6) 
Variable I YBASE (1)1 YBASE (2) 1 Ybase (3) 

Card column 10 20 IQ. 1^2 ^ 

Format type F F F F I F 1 F I F 


Item no. 12. (LSPEC*1) : 1 card (1615) 
Variable 

Card column 
Format type 


LOCO(NSELCT) 























Item no. 13. 1 card (8F10.6) 

Variable Ml Ql ) 

Card column lo" 20 \ 

Format type P p l 

Item no. 14. J cards, J as above, 8 values per card (FlO. 6 ) 

Variable [XCALB(I) JXCALBf2) XCALB(3^| 

Card column K 20 30 4J2 so 

Format type F F F F F 



Item no. 15. J cards, J as above, 8 values per car'd (8F10.6) 

Variable YCALB(l) YCALB(2) YCALbTIT 

Card column K ^ ^ 50 

Format type F F F F F 



Item no. 16. 
Variable 
Card column 

Format type 


(LSPEC-1) : 1 card (1615) 


LOCl(NSELCT) 



Item no. 17. 

Variable 
Card coliunn 
Format type 


1 card ( 8 fl 0 . 6 ) 

I M 2 1 7 


Item no. 18. 1 card (8 f10.6) 

Variable Q2(l) 1 Q2(2 

Card column W 2 

Format type F i* 


Q2 (NPARAM) 

























o 


Item no. 19 (LCHEK>“1) : J cards , J=1+INT (N/8) , 8 values per card (8F10.6) 

Variable 
Card column 
Format type 



Item no. 20. 
Variable 
Card column 
Format type 


(LCHEK=1) : J cards, J as above, 8 values per card (8P10.6) 




A. 6 DESCRIPTION OF OUTPUT 


The first output item consists of a banner page, and the 
card images of the input data, the latter only if LECHO = 1. 

The second item is a page headed by the job title, listing; 

1, the input parameters N, A, B and NPARAM relevant to the 
actual calculaton; 

2. the straining option selected and the classifications 
of the straining points selected. 

The third item is the results of the computations on the 
base solution. These include: the Mach number MO, values of 
the perturbation parameters of the base solution QO(k) , K = 1, 
NPARAM, and the critical value of the dependent parameter, the 
locations of the minimum, maximum and critical points, and the 
locations of the invariant points. These results are then 
repeated for each of the NPARAM calibration solutions. 

Results for unit straining of XBASE, and unit perturbations 
of the dependent variables are the fourth output item. This 
is done only if LUNIT = 1. 

The fifth item (repeated for each case computed) summarizes 
the results of the perturbation calculation. The Mach number, 
the values of the perturbation paramters , and the critical values 
of the variable are printed first, followed by the locations of 
the minimum, maximum, and critical points in the perturbation 
solution and comparison solution (if any) . Next, a legend is 
printed providing the maximum, minimum, and critical values 
of the dependent parameter and corresponding point symbols. 

Print symbols are also provided for the printer plot of the 
perturbation solutions P, comparison solution C (if any) , and a 
common symbol ($) when there is agreement within printer plot 
accuracy between the two. Then follows a table listing XBASE, 
YBASE, XPERT (the strained coordinate) , and YPERT (the computed 
value of the dependent variable). If LCHEK = 1, three additional 
columns list XCHEK (the computed solution) at the points given 
by XCHEK. This allows direct numerical comparison of YPERT with 
YCHEK, since the values of XPERT and XCHEK will not in general 
coincide. A printer plot is then provided of the perturbation 
result, together with the comparison solution (if any) . 

The final item is a table which summarizes the perturbation 
parameters for the base, calibration, and all predictive solu- 
tions . 
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A. 7 ERROR MESSAGES 


NUMBER OF CRITICAL POINTS IN 
BASE AND CALIBRATION SOLUTIONS 
ARE UNEQUAL - CALCULATION ENDED 

This message will be printed if critical points are speci- 
fied in straining (LSPEC = 0) and the number of critical points 
in base and calibration solutions are unequal. The remedy is 
to avoid use of critical points in straining, or to use base 
and calibration solutions having equal numbers of critical, 
points . 


NUMBER OF CRITICAL POINTS 
SELECTED EXCEEDS NUMBER 
ACTUALLY LOCATED - CALCULATION 
ENDED 

This message will be printed if more critical points are 
specified in straining (LSPEC = 0) than the number located by 
the program. The remedy is to specify a number of points less 
than or equal to the actual number. 


ORDER OF SPECIFIED POINTS IN 

BASE AND CALIBRATION SOLUTIONS 

DOES NOT CORRESPOIND - CALCULATION ENDED 

This message will be printed if the fixed points specified 
(LSPEC = 0) occur in a different sequence in the base and caliba- 
tion solutions. The remedy is to use base and calibration solu- 
tions having the same qualitative features. 


A. 8 SAMPLE CASE 

The sample case presented in this section provides some 
example results of perturbation calculations and comparisons 
with 'exact' nonlinear solutions for a multiple-shock flow for 
which partial results were provided in figure 3 of the main 
text. The calculation is for the simultaneous four-parameter 
(M^,i,t,a) perturbation of strongly-supercritical full potential 
flows past a cascade of blades having NACA four-digit profiles. 
The base flow is for oncoming Mach number M^ = 0.780, thickness 
ratio T = 0.110, gap-to-chord ratio t = 3.2, and oncoming inflow 
angle a 0.3°. The calibration flows to account for perturba- 
tions in these parameters are at the following values of these 
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parameters) (M^,Trt,a) = (0.790,0.110,3.2,0.3), (0.780,0.120,3.2, 

0. 3), (0.780,0.110,3.0,0.3), (0.780,0.110,2.3,0.5). Perturba- 

tion results have been determined for 19 different solutions, 
which are summarized in the summary table after the print output. 
Results for several of those cases are presented here. 

The input data is tabulated in figure A.l, with item numbers 
corresponding to those identified in Section A. 5.1 and A. 5. 2. 

The first card, item 1, provides the title of the run. The 
next card, item 2, indicates that the input deck will not be 
printed (LECHO = 0) . The next card, item 3, indicates that 
there are 191 points (N = 191) at which data v/ill be input for 
the base, calibration, and comparison solutions; that there 
will be 19 cases (NCASE = 19) for which perturbation solutions 
are to be computed; that the number of parameters to be 
perturbed are 4 (NPARAM = 4) ; that there will be three invariant 
points (NSELCT =3) in addition to the end points; that the 
invariant points will be located by the program (LSPEC = 0) ; 
that the information regarding the unit perturbation will be 
printed (LUNIT = 1) ; that there will be a comparison of the 
perturbation results with the exact solution (LCHEK = 1) ; and 
that there will be plots by a peripheral device of the output 
(LPLOT = 1) . The next card, item 4, specifies that the three 
invariant points to be selected by the program are to be: the 
medium point (LSELECT(l) = 2) , i.e., the stagnation point; the 
first critical point (LSELCT(2) = 3) , i.e., the 1st shock point 
found when moving forward on the bottom surface from the 
leading edge; and the second critical point (LSELCT(3) = 4) , 

1. e., the end shock point. 

The next card, item 5, indicates that the dependent vari- 
able for print output will be symbolized by a 'CP' denoting 
the pressure coefficient. The next card, item 6, indicates 
that the parameters to be varied are "PARNAM(l) = MACH NO., 
PARNAM(L) = TAU, PARNAM(3) = PITCH, PARNAM ( 4 ) = ALPHA 1 . The 

next card, item 7, indicates that the coordinates of the data 
points to be read in will start at x = 1.0 on the lower surface 
and will end at x = 1.0 on the uppper surface, i.e., A = -1, 

B = 1. The next card, item 8, provides the oncoming Mach number 
of the base flow MO = 0.780. The next card, item 9, provides 
the base flow values of the parameters to be perturbed: Q0(1) = 
0.780, Q0(2) = 0.110, Q0(3) = 3.2, Q0(4) = 0.30. The following 
25 cards, item 10, provide the 191 base flow values of the 
surface coordinates XBASE (I) , I = 1, N, while the next 5 cards, 
item 10, provide the 191 base flow values of the surface coordi- 
nates XBASE (I) , I = 1, N, while the next 25 cards, item 11 
provide the 191 base flow values of the dependent variable 
(pressure coefficient) YBASE(I) , I = 1, N. The next card, item 
13 provides the values of the oncoming Mach number Ml and the 
value of the Ist-perturbation parameter in the first calibration 
solution. Items 14 and 15, which correspond to the arrays 
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of 191 points of coordinates XCALB(I) , I = 1, N and dependent 
variable YCALB(I), I = 1, N, are provided as for the base flow 
in two sets of 25 cards each. Items 13, 14, and 15 are then 
repeated 3 more times corresponding to the total of NPARAM = 4 
calibration solutions required. Items 17, 18, 19 and 20, of 
which there are 19 sets corresponding to the 19 cases to be 
studied, provide analogous information as items 8, 9, 10, and 
11 of the base flow, but now refer to the 'exact' nonlinear 
results. Items 19 and 20, which correspond to the XCHEK(I), 

I = 1, N and YCHEK(I), I = 1, N arrays, respectively, have 
of course been previously compted at the indicated values of 
the perturbed parameters, and are included here for comparative 
purposes to enable assessment of the perturbation results. 

Figure A, 2 provides an abbreviated print output for 
sample case, while figure A. 3 provides an abbreviated plot 
output of the results for the 19 cases, and display the base 

( ), calibration (••••) r perturbation (****), and 'exact" 

nonlinear ( ) flow solutions. In these plots, the calibra- 

tion solutions denoted "calibration no. 1 of 4", etc. corre- 
spond to those indicated in the summary table provided at the 
end of the print output for this case shown in Figure A. 2. 
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TfcST CASE - 4 parameter PERTURBATION Of A SUPERCRITICAL CASCADE FLOW 
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Figure A. 2- Abbreviated print output for sample case 




• test case * 


PROORAM PERTURB 

CALCULATES NONLINEAR MULTIPLE-PARAMETER 
CONTINUOUS OR DISCONTINUOUS 
PERTURBATION SOLUTIONS 
WHICH REPRESENT CHANGES IN EITHER 
GEOMETRY OR FLOW CONDITIONS 
BY EMPLOYING A STRAINED-COORDINATE PROCEDURE 
UTILIZING UNIT PERTURBATIONS DETERMINED FROM 
PREVIOUSLY CALCULATED 
••BASE" AND "CALIBRATION" SOLUTIONS 
DISPLACED FROM ONE ANOTHER BY SOME REASONABLE 
CHANGE IN GEOMETRY OR FLOW CONDITION 

WRITTEN RY 

JAMES P. ELLIOTT AND STEPHEN S. STAHARA 

NIELSEN engineering AND RESEARCH* INC. 
MOUNTAIN VIEW* CALIFORNIA 



Figure A. 2- Continued 


TtST CASE - 4 PARAHETtH PERTURBATION OF A SUPERCRITICAL CASCADE FLOW 


LIST OF INPUT parameters 
N = 191 

A s *1.0 8 = 1.0 

NPARAM « 4 


straining OPTIONS 
NUMBER OF FIXED POINTS: S 

FIXED POINTS WILL BE AUTOMATICALLY DETERMINED 
BY THE PROGRAM FOR ALL SOLUTIONS AS FOLLOWS: 

TWO END POINTS 
POINT OF MAXIMUM CP 
CPCRIT (1ST POINT) 

CPCRIT (2ND POINT) 


Ov 

tn 



Figure A. 2- Continued 


WESULTS Of computations ON BASE SOLUTION: 
MACH NUMBER* 

VALUES OF PERTURBATION PARAMETERS* 
CRITICAL VALUE OF CP: 


II 

o 

X 

7800 


00(1) = 

.7800 

(MACH NO. 

00(2) = 

.1100 

( TAU 

00(3) > 

3.2000 

( PITCH 

00(4) = 

.3000 

( ALPHAI 

CPCRIT 

s -.4940 



...LOCATIONS OF MIN,* MAX.* AND CRITICAL PTS 
(* DENOTES POINT ON LOWER SURFACE) 


MINIMUM AT X = .A220* 

MAXIMUM AT X := 0.0000 
2 CRITICAL POINT (S) : 
1ST AT X = .5001* 

2ND AT X = .4273 


(POINT NO. S3) 
(POINT NO, 96) 

(after POINT NO. 
(AFTER POINT NO. 


...LOCATION OF FIXED POINTS 
(* DENOTES POINT ON LOWER SURFACE) 

XFIX(l) = 1.0000* 

XFIX(2> » .5001* 

XFIX(3) = 0.0000 
XFIX(4) 3 .4273 

XFIX(S) > 1.0000 


46) 

139) 



Figure A. 2- Continued 


HtSULTS OF COMPUTATIONS ON 1ST CALIBRATION SOLN; 


I....MACH NUMBERf 

VALUES OF PERTURBATION PARAMETERS, 
CRITICAL VALUE OF CP: 

Ml = .7900 

DENOTES PERTURBATION FROM BASE VALUE) 
*•01(1) « .7900 (MACH NO.) 

01(2) s .1100 ( TAU ) 

QK3) * 3.2000 ( PITCH ) 

01(4) = .3000 ( ALPHAl ) 

CPCRIT = -.4638 


...LOCATIONS OF MIN., MAX,, AND CRITICAL PTS 
(* DENOTES POINT ON LOWER SURFACE) 


MINIMUM AT X = .5197* 

MAXIMUM AT X = 0.0000 
2 CRITICAL POINT(S) : 
1ST AT X = ,5788« 

2ND AT X = .5149 


(POINT NO, 47) 
(POINT NO. 96) 

(AFTER POINT NO. 
(AFTER POINT NO. 


...LOCATION OF FIXED POINTS 
(* DENOTES POINT ON LOWER SURFACE) 


XFIX(l) 

XFIX(2) 

XF1X(3) 

XFIX(4) 

XFIX(5) 


1 . 0000 * 

.5788* 

0.0000 

.5149 

l.OOOO 


43) 

144) 


ON 
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00 WESULTS OF COMPUTATIONS ON 2ND CALIBRATION SOLN: 


....MACH NUMBER « 

VALUES OF PERTURBATION PARAMETERS. 
CRITICAL VALUE OF CP: 

Ml = .7800 

(** DENOTES PERTURBATION FROM BASE VALUEJ 


OKI) = 

.7800 

(MACH NO. 

**Ql(2) X 

.1200 

t TAU 

01(3) = 

3.2000 

( PITCH 

01(4) = 

.3000 

( ALPHAI 

CPCRIT = 

-.4940 



...LOCATIONS OF MIN., MAX., AND CRITICAL PTS. 
<* DENOTES POINT ON LOWER SURFACE) 


MINIMUM AT X = .5197* 

maximum at X - 0.0000 
2 CRITICAL POINT(S) : 
1ST AT X = .5774* 

2ND AT X X .5175 


(POINT NO. 47) 

(POINT NO. 96) 

(AFTER POINT NO. 43) 
(AFTER POINT NO. 144) 


...LOCATION OF FIXED POINTS 
(• DENOTES POINT ON LOWER SURFACE) 


XFlX(l) X 1.0000* 
XFIX(2) X .5774* 
XFIXO) X 0.0000 
XFIX(4) X .5175 
XFIX(5) > I. 0000 
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HtSULTS OF COMPUTATIONS ON 3 RD CALIBRATION SOLN 


....MACH number. 

VALUES OF PERTURBATION PARAMETERS. 
CRITICAL VALUE OF CP: 

Ml = .7800 

(*• denotes PERTURBATION FROM BASE VALUE) 


Old) = 

.7800 

(MACH NO. 

QK2) = 

.1100 

( TAU 

**(J1(3) = 

3.0000 

( PITCH 

Ql(4) = 

.3000 

( ALPHA 1 

CPCRIT 3 

-.4940 



...LOCATIONS OF MIN.. MAX.. AND CRITICAL PTS 
<* DENOTES POINT ON LOWER SURFACE) 


MINIMUM AT X = .4382* 

MAXIMUM AT X s 0.0000 
2 CRITICAL point (S) ! 
1ST AT X = .5165* 

2ND AT X s .4520 


(POINT NO. 52) 
(POINT NO. 96) 

(AFTER POINT NO. 
(AFTER POINT NO. 


...LOCATION OF FIXED POINTS 
(* DENOTES POINT ON LOWER SURFACE) 

XFIX(l) 3 1.0000* 

XFIX(2) » .5165* 

XFIX(3) s 0.0000 
XFIX(4) 3 .4520 

XFIX(5) = 1.0000 


o\ 

VO 


47) 

140) 
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RESULTS or COMPUTATIONS ON 4TH CALIBRATION SOLN ! 


MACH NUMBER t 

values of perturbation PARAMETERSf 
CRITICAL VALUE OF CP: 


Ml = 

.7800 


DENOTES PERTURBATION 

FROM BASE VALUE: 

Q1 (1) 

a .7800 

(MACH NO.) 

01 (2) 

= .1100 

( TAU ) 

01(3) 

= 3.2000 

( PITCH ) 

**Ql (4) 

a ,5000 

( ALPHA 1 ) 


CPCRIT = -,^940 


LOCATIONS OF MIN., MAX.t AND CRITICAL PTS. 

(* DENOTES POINT ON LOWER SURFACE) 

MINIMUM AT X = ,4058* (POINT NO. 54) 

MAXIMUM AT X a 0,0000 (POINT NO. 96) 

2 CRITICAL POINT(S) : 

1ST AT X a ,4903* (AFTER POINT NO, 48) 

2ND AT X « ,4360 (AFTER POINT NO. 139) 


,,. LOCATION OF FIXED POINTS 
(* DENOTES POINT ON LOWER SURFACE) 

XFIX(l) a 1,0000* 

XFIX(2) a ,4903* 

XFIX(3) a 0,0000 
XFIX(4) a ,4360 
XFIX(5) a 1.0000 



UNIT PEWTUHBATION OF CP AND UNIT STRAINING OF XBASE 
for calibration solutions I THROUGH 4 
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iHfl .9955 .9962 -.4223 .9962 2.3418 .9957 .0160 .9956 .0261 

li^9 .9973 .9978 -1.4918 .9978 1.1609 .9975 .0313 .9974 .0283 

190 .9987 .9989 -7.6871 .9989 -.9395 .9987 .0648 .9987 .0275 

191 .9995 .9996 -24.6958 .9996 -1.9799 .9995 .1357 .9995 -.0299 
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• OUTt»UT FOR CASE NO. 1 OF 19 * 


.MACH number. 

values of perturbation parameters, 
critical value of CP: 

M2 « .7750 

Q2<ll = .7750 (MACH NO.) 

02(2) = .1200 ( TAU ) 

02(3) s 3.0000 ( PITCH > 

02(4) s .5000 ( ALPHAl ) 


CPCRIT = -.5095 

H* 

.....LOCATIONS OF HIN.. MAX,. AND CRITICAL PTS. 
^ (* denotes POINT ON LOWER SURFACE) 

(D 

^ PERTURBATION SOLN: 


to 

I 

n 

o 

3 

ft 

H- 

3 

3 

(D 

pj 


MINIMUM AT X = ,4771* (POINT NO. 52) 

MAXIMUM AT X = 0.0000 (POINT NO. 96) 

2 CRITICAL POINT(S) : 

1ST at X a ,5420* (AFTER POINT NO. 48) 

2ND AT X a ,5065 (AFTER POINT NO. 139) 

COMPARISON SOLN: 

MINIMUM AT X « .4870* (POINT NO. 49) 

MAXIMUM AT X B 0.0000 (POINT NO. 96) 

2 CRITICAL POINT(S) ( 

1ST AT X « .5493* (AFTER POINT NO. 45) 

2ND AT X s .5066 (AFTER POINT NO. 144) 


FINAL PRINTOUT AND GRAPHICAL OISPLAT OF CP 


PT 

XBASE 

CPBASE 

XPERT 

CPPERT 

xchek 

CPCHEK 

CPPINT H- 

1 

.9945 

.6608 

.9995 

.8210 

.9995 

.7155 

.8108 

2 

.9987 

.5751 

.9988 

.6594 

.9987 

.6231 

.6506 

3 

.9974 

.5086 

,9976 

.5706 

.9974 

.5555 

.5626 

4 

.9956 

.4543 

.9960 

.5145 

.9956 

.4972 

.5050 

5 

.9932 

.4077 

.9938 

.4637 

.9932 

.4469 

.4533 

6 

.9903 

,3659 

.9912 

.4197 

.9903 

.4022 

.4081 


H a MAXIMUM VALUE OF CP * 1.1517 

L a MINIMUM VALUE OF CP « -1.0527 

* = CRITICAL value of CP « -.5095 

P a VALUE OF CP PREDICTED BT PERTURBATION SOLUTION 

C a VALUE OF CP IN COMPARISON SOLUTION 
« a agreement BETWEEN P AND C 

• L 

PC * 

PC ♦ 

PC * 

PC • 

s * 

PC * 



S>P-S9 .3274 .98ftl .3792 .9869 .3619 .3667 
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Figure A. 2- Continued 


* (tUlPUT POP CASE NO. 1<J OF 19 • 


.MACH NUMBER. 

VALUES OF PERTURBATION 
CRITICAL VALUE OF CP: 

parameters 

M2 = .7900 


02(1) = .7900 

(MACH NO. 

0212) = .1200 

( TAU 

02(3) = 3.0000 

( PITCH 

02(A) = .2000 

( ALPHAl 

CPCRIT = -.4636 



...LOCATIONS OF MIN.. MAX.. ANO CRITICAL PTS. 
(• DENOTES POINT ON LOmER SURFACE) 


PERTUReATION SOLN: 


MINIMUM AT 

X = 

.6154* 

(POINT 

NO. 

51) 

MAXIMUM AT 

X = 

.0000* 

(POINT 

NO. 

96) 

2 CRITICAL 

POINT(S) : 




1ST AT 

X = 

.6726* 

(AFTER 

POINT 

NO. 48) 

2ND AT 

X = 

.6166 

(AFTER 

POINT 

NO. 138) 

COMPARISON SOLN: 






MINIMUM AT 

X = 

.6939* 

(POINT 

NO. 

36) 

MAXIMUM AT 

X = 

0.0000 

(POINT 

NO. 

96) 

2 critical 

POINT (S) : 




1ST AT 

X = 

.7230* 

(AFTER 

POINT 

NO. 34) 

2ND AT 

X = 

.6903 

(AFTER 

POINT 

NO. 155) 


FINAL PRINTOUT AND GRAPHICAL DISPLAY OF CP 


H 

s MAXIMUM 

VALUE 

OF 

CP = 

1.1565 

L 

s MINIMUM 

VALUE 

OF 

CP = 

-1.1126 

• 

= CRITICAL 

VALUE 

OF 

CP = 

-.4638 


P = VALUE OF CP PREDICTED BY PERTURBATION SOLUTION 
C = VALUE OF CP IN COMPARISON SOLUTION 
t = AGREEMENT BETWEEN P AND C 


KD 


PT 

XBASE 

CPBA5E 

XPERT 

CPPERT 

XCHEK 

CPCHEK 

CPPINT H. 

1 

.9995 

.6608 

.9997 

,7696 

,9995 

.7942 

.7631 

? 

.9987 

.5751 

.9992 

.7513 

.9987 

.7350 

.7120 

.3 

.9974 

.5086 

.9963 

.6791 

.9974 

.6519 

.6249 

A 

.9956 

.45<»3 

.9971 

.6098 

.9956 

.5605 

.5545 

S 

.9932 

.4077 

.9956 

.5560 

.9932 

.5367 

.4993 


.9903 

.3659 

.9938 

.5100 

.9903 

.5007 

.4500 


— * L 

CP * 

$ * 

CP * 

4 ♦ 

CP * 

CP * 
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c 


PRO0B*N PERTRB MNPUTtOUTPUTtTAPClf T«PCS<INPUT,TAPE6>0UTPUT« 
% T»PEJ5» 


MAIN PROGRAM PERTURB 


CALCULATES CONTINUOUS OR DISCONTINUOUS NONLINEAR PERTURBATION 
solutions NHICH represent a multiple-parameter CHANGE IN EITHER 
GEOMETRY OR FLOW CONDITIONS SV EMPLOYING A STRAINEO-COORDINATE 
PROCEDURE, the METHOD UTILIZES UNIT PERTURBATIONS. DEIERHInEO 
For each parameter from a previously calculated common base 
SOLUTION, and a CALIBRATION SOLUTION DISPLACED FROM IT BY SOME 

Reasonable change in the relevant parameter, to predict new 

NONLINEAR SOLUTIONS OVER A RANGE OF PARAMETER VARIATION. 

THIS CURRENT VERSION OF THE PROCEDURE IS CONFIGURED TO PREDICT AND 
PLOT pressure COEFFICIENTS ON A BLADE OR AIRFOIL SURFACE. AND CAN 
ACCOUNT FOB THE MOTION OF 

<I> ONE OR MORE CRITICAL POINTS ISHOCK POINTS). 

<21 A stagnation POINT. 

<3T A MAX IMUM-SUCT ION-PRESSURE POINT. 

OR SIMULTANEOUSLY FOR ANY COMBINATION OF THESE. 

the PROORAM IS ALSO CONFIGURED TO COMPARE THE PERTURBATION- 
PREDICTED SOLUTIONS WITH THE CORRESPONDING EXACT SOLUTIONS 
OBTAINED BY EMPLOYING THE SAME EXPENSIVE COMPUTATIONAL 
PROCEDURE USED TO DETERMINE THE BASE AND CALIBRATION SOLUTIONS. 

SEE the SUBROUTINE INPUT FOR DETAILS. 

N ■ NO. OF POINTS IN SURFACE PRESSURE DISTRIBUTION - ASSUMED EQUAL 
FOR BASE. CALIBRATION. AND PREDICTED DISTRIBUTIONS. 

NOTE) N .LE. 200. 

NPARAH > NO. OF PARAMETERS VARIED. NOTE) NPARAM .LE. B. 


BASE 

MACH NO. MO 

PERT. PARAHS. 00 ( 1 ) .00<2) • . . . 00(11 


KTH calibration 


.Q0<K-|).Q|,Q0IK*11 < 


H2 ■ ONCOHINO HACH NO. OF PREDICTED FLOW. 

02(1) .0212) I... ■ VALUES OF PERTURBED PARAMETERS IN PREDICTED 
FLOW. 

COORDINATE STRAINING IS PIECEWISE LINEAR WITH END POINTS AND ONE 
OR HORE USER-SELECTED INTERIOR POINTS HELD INVARIANT. 

THE PROGRAM LOCATES MINIMUM, MAXIMUM, AND ALL CRITICAL POINTS 
ISHOCK POINTS) IN THE BASE AND CALIBRATION SOLUTIONS. AND STORES 
these in the arrays XLOCO AND XLOCl (IT IS ASSUMED THAI THE NUMBER 
OF critical points does not exceed FOUR) XS FOLLOWS 

BASE CALIRRATION 


XLOCOd) T MINIMUM PT. XLOCHI) > MINIMUM PT, 

XLOCO(2) > MAXIMUM PT. XL0CK2) s MAXIMUM PT. 

XLOCOI3) * CPITICAL PT. NO. I XLOCIO) * CRITICAL PT. NO. 1 


XLOCOI&) ■ CRITICAL PT. NO. A 




c 

IS tHE TRCC STIICAM MCH NUNSERt AND IGRAO C0RRESP0M>S TO POSITIVE 

MAINIIS 

c 

PRESSURE ORADIENT <*II, 

MAINI16 

c 


HAINIIT 


VCRlT(Zl>2.0*l<(2>0*0.*«2*«2l/2.A)**(I.A/O.Al-|.OI/(l.**Z*»2) 

MAINIIO 


IGRA0>) 

HAINII9 

c 


MAINI20 



MAIN121 

c 


MA1N122 

c« • • • 

.PRIRT RANNER PAGE. 

MAIN123 

c 


MAINI2A 


call BARNER 

HA1M12S 

c 


MAIN126 

c..«. 

•READ LECHO AND ECHO INPUT DECK IE LECHO .EO. I. 

MAINI27 

c 


HAIN12B 


CALL ECHINP 

HA1N129 

c 


MAIN130 

c« • « • 

•INPUT CONTROL* GEOMETRY* AND STRAINING PARAMETERS. 

MA1N13I 

c 


MAIN! 32 


CALL INPUT Ul 

MAINI33 

c 


HA1N13A 

C« f » • 

•MRITE TITLE and INPUT PARAMETERS. 

HA1N13S 

c 


MAIN! 36 


WRITE 16.10001 TITLE 

MA1NI3T 


WRITE I6.10I0> N*A.B*NPARAM 

MAIN! 38 

c 


MA1N139 


NEIA«NSELCT*2 

MAINUO 


NSEG«NFIX-I 

MAINlAl 

c 


MAINU2 

c.... 

.PRINT INFORMATION REGARDING STRAINING TO BE USED. 

HAINIA3 

c 


MAINlAA 


WRITE <6*1020) NFIX 

MAINlAS 


IF ILSPEC .EG. 0) GO TO 10 

MA1NI46 


WRITE 16*10301 

MAINIAT 


GO TO 30 

HAIN1A6 

10 

continue 

MA1N1A9 


WRITE (6*I0A0I 

MAINlSO 


00 20 1>I*NSELcT 

MAINISI 


IF ILSELCTlI) .EG. 1> WRITE <6.I0S0) VRAM 

MAIN1S2 


IF ILSELCT(I) .EG. 2) WRITE J6.1060) VRAM 

HAIN1S3 


IF ILSELCTd) .LE. 2) 00 TO 20 

MAINISA 


LC0RR«I 

MAINISS 


LPR«LSELCT 1 1) -2 

MAIRlSb 


WRITE <6*1070) VNAM*ORD<LPR) 

MAINIST 

20 

CONTINUE 

HA1N158 

30 

CONTINUE 

MAIR1S9 

C 


HA1NI60 

c.... 

•BEGIN CALCULATIONS ON BASE SOLUTION. 

MAINlA) 

c 


HAIN162 


CALL INPUT <2> 

MAIN163 


rCRO>VCRIT|M0) 

HAIN16A 


WRITE 16.1060) HEADO 

MA1N16S 


WRITE 16.1090) VNAR 

MAIN166 


WRITE (6.1100) HO 

MAIN167 


IF (NPAHAM .EQ. 1) WRITE (6*1110) 00 ( 1 ) .PARNAH | p 

MA]Nl68 


IF (NRARAM .6T. 1) WRITE (6.1120) (K,(}0 (K) .PAHNAM (K ) .K = 1 .NPAPAM) 

MAINI69 


WRITE (6.1130) VNAM.YCPO 

haEn] 70 

c 


maim 71 


c» 

to 


c. 

c 

c 


c 

c. 

c 

c 


c 

c. 

c 


c 

c 

c 

c 

c 

c, 

c 


...NORMALIZE X coordinates ANO LOCATE MINIMUM. MAXIMUM. AND CRITICAL 

MA 

INIT2 

POINTS FOR BASE SOLUTION. 



HA 

INI73 




MA 

INIT4 

CALL scale (N.XBASE.I.A.B) 



HA 

INITS 

CALL LOCATE (N. XBASE .VBASE.YCRO. IGRAD.LMNO.LNXO* 

ncro.lcro 

•XLOCO) 

HA 

IN176 

YeCMlNiYBASE(LMNO) 



MA 

1N177 

YBCHAXxYBASE (LMXQ l 



HA 

IN176 

WRITE (6.1140) 



MA 

IN1T9 

WRITE (6.11S0) 



HA 

INIBO 

call UPLOW (A*B.XLOCO*6*NCRO*2*XDUT*FLAG) 



HA 

INIM 

WRITE (6.1160) XOUTID.FLAGtl) •L(«RO*XOUT<Z) *FLAG(2) .LMXO 


HA 

1N1M2 

IF (NCRO .GT. 0) WRITE (6. 1 170) NCRO. 



MA 

IN1M3 

6 (OROd ) .XOUT (1*2) .FLAG(I«2) .LCRO ( 1 ) * I >1 .NCRO) 


HA 

INIn4 




MA 

INlnS 

...LOAD SELECTED STRAINING POINTS INTO FIXED-POINT 

array for 

BASE 

ha 

IN1(I6 

SOLUTION. 



HA 

1N1B7 




HA 

IN1B6 

XFIX0<I>-0.0 



ma 

1N1B9 

XFIX0(NFIX)b1.0 



HA 

INI90 

IF (LSPEC .eg. I) WRITE (6.I1B0) 



MA 

IN191 

DO 50 I-I.NSELCT 



MA 

INI42 

IF (LSPEC .EG. 0) GO TO 40 



MA 

IN193 

XF I XO ( I * 1 ) xXBASE (LOCO ( I ) ) 



HA 

IN194 

WRITE (6. 1190) LOCO(I) 



HA 

1N195 

GO TO SO 



MA 

INI96 

40 CONTINUE 



H* 

INI97 

XFIX0(I*1)>XL0C0<LSELCT<1) I 



HA 

IN 196 

50 CONTINUE 



Ha 

IN 199 




MA 

1N200 

...arrange selected fixed POINTS IN A MONOTONE SEOUENCE. 


HA 

IN201 




MA 

IN202 

CALL SORT (NFIX.XFIXO.ISEOO) 



ha 

IN203 

WRITE (6.1200) 



HA 

IN20A 

WRITE (6.11S0I 



HA 

1N205 

CALL UPLOW (A.B.XFIXO.BiNFIXiXOUT.FLAG) 



ma 

IN206 

WRITE (6.1210) II.XOUT(l).FLAG<I)*I>|.NFlX) 



MA 

IN207 




HA 

(N2nB 

...END CALCULATIONS OH BASE SOLUTION. 



MA 

1N209 




MA 

IN210 




ha 

1N21I 




Ha 

IN212 

...BEGIN CALCULATIONS OR CALIBRATION SOLUTIONS. 



ha 

1N213 




HA 

(N214 

DO 160 Kel.NPARAM 



ma 

IN21S 

CALL INPUT (3) 



HA 

IN2I6 

M1SAVE(K)»M1 



MA 

IN217 

01SAVE(K)>01 



ha 

[N2IB 

TCR1*YCRIT(HI) 



ha 

[N2I9 

DELI (K)>01-00IK) 



ma 

IN220 

call copy M.N.K.XCALB.XCSAVE) 



MA 

IN221 

CALL COPY (l.N.K.YCALB.YCSAVE) 



HA 

IN222 

IF (NPaRAM .eg. 1) WRITE (6. 1060) HEAOl 



HA 

IN223 

IF (NPARAH .gT. 1) WHITE (6.1220) 0RD(K).HEAD1 



HA 

[N224 

WRITE (6.1090) VNAM 



HA 

IN225 

WRITE (6.1230) Ml 



HA 

IN226 

IF (NPaRAM .GT. 1) WRITE (6.1240) 



ha 

IN22T 

no 6D KKsUNPARAM 



HA 

IN22« 


VO 

o 


c 

c 

c 

c 


c 

c 

c 


c 

c. 

c 

c 

c 

c 

c 

c 

c 

c 

c. 

c 

c 


IF INPARAM ,E0. 1) WRITE (Ail?50l OliPARNAM(l) 

IF (NPARAM .GT. J .AND. K« .tO. A) WRITE (6.1ZA0) KK.Ql .PARMAM (KK ) 
IF (KK .NE. K1 WHITE (<>ilZT01 KA ,00 ( KK ) , PARNAM | KK I 
60 CONTINUE 

WRITE 16,11301 VNAHfYCRl 

...norhaliee X coordinates and locate minimum, maximum, and critical 

POINTS FOR KTH CALIRRATION SOLUTION. 

CALL SCALE IN,XCALBtl«A,B> 

CALL LOCATE (N.XCALB, TCALBi TCP 1 , IGRAO,LHN| ,LMX 1 ,NCR1 ,LCR I • XLOCl I 
YTM1N>YCALB<LMN1) 

VTMAX*YCALR<LMk1) 

IF (TTMIN ,LT, YBCMINI YBCMIN-YTMIN 
IF lYTMAX ,GT, YBCMAXI TBCMAX*YTMAX 
WRITE (6,1140) 

WRITE (6,11S0) 

CALL UPLOW )A,B,XL0CI«6,NCR1*2,X0UT,FLAG) 

WRITE (6,1160) X0UT(1).FLA0(1),LMNI,X0UT(21,FLA6(Z),LMX1 
IF (NCR! ,GT, 0) WRITE (6,HT0) NCRt, 

% (ORO(t) ,X0UT(I*2) , FLAG (1*2) ,LCR1 ( I ) , I ■ 1 ,NCR1 ) 


•CMECA for invalid straining SPECIFICATION If LSPFC * 0. 


IF (LSPEC ,EO. 1) GO TO BO 
I COUNT. 0 

DO TO 1.I»NSELCT 
IF (LSELCT(I) ,LE. 2) 00 TO TO 
I COUNT. I COUNT *1 
IF (NCRO ,NE. NCRl) LTERM.I 
70 CONTINUE 


,ST0P EXECUTION IF CRITICAL POINTS ARE TO BE USED IN STRAINING AND 
number of critical POINIS IN BASE AND CALIBRATION SOLUTIONS ARE 
UNEQUAL, 


IF (LTERM .EQ. 1> GO TO 900 


.STOP EXECUTION IF NUMBER OF CRITICAL POINTS SELECTED EXCEEDS 
NUMBER actually LOCATED. 


IF (ICOUNT ,GT, NCRO) GO TO 90S 


80 CONTINUE 


.LOAD SELECTED STRAINING POINTS INTO FIXED-POINT ARRAY FOR KTH 
calibration SOLUTION. 


XFIXI (I).O.O 
XFIXl (NFIX) >1.0 

IF (LSPEC ,EO. 1) WRITE (6,11B0) 
00 100 I»1,NSELCT 
IF (LSPEC .EO. 01 GO TO 90 
XFIXI (I*1IcXCAlB(LOCI(I)I 
WRITE (6,1190) LOCIII) 

GO TO 100 
90 continue 


XFIXI (I.I).XLOCI ILSELCT(I) I MAINZR6 

100 continue MAIN2A7 

^ MAIN2BB 

ARRANGE SELECTED FIXED POINTS IN A MONOTONE SEQUENCE. MAIN2H9 

C MA1N290 

call sort (NFIX, XFIXI, ISEOI) MA1N291 

WRITE (6,1?00) MA1N292 

WRITE (6.1IS0) HAIN293 

CALL UPLOW (A.B,XFIX1,8«NFIX,X0UT,FLAG) MA1N294 

WRITE (6,1210) I1,X0UT(1),FLA6(1),1.|,NFIX) MA1N29S 

C MA1N296 

C STOP EXECUTION IF ORDER OF OCCURRENCE OF CRITICAL POINTS IN BASE MAIN29T 

C AND calibration SOLUTIONS DOES NOT CORRESPOND. HAIN298 

C MA1N299 

IF (LSPEC .EO. I) GO TO 120 MAIN300 

DO no I>1,NF1X MAIN301 

IF (ISEOO(I) ,NE. ISEOKII) GO TO 910 MAIN302 

110 continue Mx1N3n3 

120 continue MXIN304 

C MA1N30S 

C COMPUTE COEFFICIENTS IN KTH UNIT STRAINING OF XBASE MAIN306 

C MAIN307 

C XSTR > C(K.l) . 0(K.t)*XBASE. I>1.2, ... .NSEG. MAIN30B 

C MAIN309 

C where NSEG IS THE NUMBER OF LINEAR SEGMENTS. MAIN310 

C M41N311 

DO 130 I'ltHSEG MA1N3I2 

CNUM.XFIX1 (I)*XFIX0(1»D-XF1X1 H»1)»XFIX0(I) MA1N313 

DNUM.XFlXI(I«))-XFIX|(II MA1N314 

OENOMexFIXO(I‘l)-XFI40(l) M«|N 3|5 

C(K,I)=CNUM/DEN0M MAIN3I6 

0(K,l).0NUH/DEN0M MAIN3I7 

130 continue MAINJIB 

C MA1N3I9 ' 

C determine KTH UNIT STRAINING OF XBASE. MAIN320 

C MAIN321 

call strain (N.K, NSEG, XFIXO, XBASE, l.O.OELX) MAIN322 

DO 140 I«l.N MAIN323 

140 XUNIT(l). XBASE (1)»0ELX(I) MAIN324 

C HAIN325 

C interpolate CALIBRATION SOLUTION TO BASE FLOW POINTS CORRESPONDING MAIN326 

C TO UNIT straining, HA1NJ27 

c MAIN3?a 

call INTERP IN,XCALB,YCALB,XUNIT,Y1NTP) HAIN329 

C MA1N330 

C CORRECT VALUES ON EITHER SIDE OF CRITICAL POINTS. IF THESE ARE HAIN331 

C USED IN straining. HA1N332 

C MAIN333 

IF (LCORR .EO. 0) GO TO 160 MAIN334 

00 150 I.l.NCRl HAIN335 

YINTP(LCR0(I) I.YCALBILCRKin MAIN336 

YINTP(LCRO(I)«1).yCALB(LCR1 (I) M) MAIN337 

ISO continue HAIN33K 

160 CONTINUE HAIN339 

C MAIN340 

C DETERMINE THE KTH UNIT PERTURBATION. MAIN3*1 

C MAIN342 





00 170 Ia|,N 

HA 

IN343 

c 


HAIN400 


I7t> 

YUN1T|K>II«IYINTPI|)-V8«SEII))/0ELMKI 

HA 

IN34* 

C.... 

.initialize strained coordinate and perturbation solution. 

MAIN401 

c 



HA 

IN345 

C 


MAIN402 

c 


.SAVE UNIT STRAINING IP REQUIRED FOR LATER PRINTOUT. 

HA 

IN346 


do 250 lal.N 

HAIN403 

c 



HA 

IN34 7 


XPERT II )*XBASEII) 

HAJN404 



IF (LUNIT .EO. 0) GO TO ISO 

MA 

IN348 

250 

YPERTIDaYBASCIII 

MAIN405 



CALL SCALE (N.XUNIT.2< A«BI 

ma 

IN349 

C 


HAIN4n6 



Call copy iun.k.xunit.xusavE) 

HA 

1N350 

c.... 

.ADO IN contributions FROM ALL PERTURBATIONS. 

HAIN407 


laa 

CONTINUE 

HA 

IN3SI 

c 


MAIN40B 

c 



MA 

1N3S2 


DO 270 Kal.NPARAM 

HAIN409 

c 

«*• 

>ENO calculations ON CALIBRATION SOLUTIONS. 

HA 

1N3S3 


02SAVElK.ICASE)a02IK) 

HAIN4I0 

c 



HA 

IN3S4 


OEL2aQ2<K)-OOIK) 

HA1N411 

c< 

••• 


HA 

IN355 


OEL21aDEL2/DELI IK) 

HA1N412 

c 



ha 

IN3S6 


call strain in. K.NSEG.xFIXO. XBASE. 0EL21.DELX) 

MAIN413 

c. 

• •• 

PRINT UNIT PCRTURBATIONISI AND UNIT STRAININGIS) IF LUNIT ,NE. 0. 

MA 

IN3S7 


00 260 lal.N 

HAIN414 

c 



ha 

1N3SB 


XPERTIDaXPERTID.OELXII) 

HA1N4IS 



IF (LUNIT .EO. 01 60 TO 240 

HA 

IN3S9 

260 

YPERT 1 1 ) a yPERT 1 1 ) »0EL2aYUNIT IK. 1 ) 

MA1N414 



call SCALE IN.XBASE«2«A.B) 

HA 

IN360 

270 

CONTINUE 

HAIN4I7 



IRPTan 

HA 

IN361 

C 


HAIN418 



IF INPARAH .GT. 4) IRPTsl 

HA 

IN362 

C.... 

.ADJUST values near CRITICAL POINT FOR MONOTONE BEHAVIOR. 

HAIN419 



KSTARTal 

ha 

1N363 

c 


HAIN429 



KSTOPa* 

ha 

IN364 


IF ILCORR .EO. 11 CALL MONO INCRO.LCRO. XPERT. YPERT) 

MAIN421 



IF IKSTOP .GT. NPARAMI KSTOPaNPARAH 

MA 

IN365 

c 


HA1N422 



60 TO 200 

ma 

IN366 

c.... 

.LOCATE HININUH. MAXIMUM. AND CRITICAL POINTS IN PERTURBATION 

HAIN423 


190 

KSTARTaS 

HA 

IN367 

c 

SOLUTION. 

MAIN* 24 



KSTOPaNPARAH 

HA 

1N36B 

c 


HAIN425 


200 

CONTINUE 

HA 

IN369 


call scale IN.XPERT.2.A.B) 

MAIN426 



HRITE (6. 1280) VNAH 

HA 

IN3T0 


CALL SCALE IN.XPERT. l.A.B) 

HAtNA27 



IF INPARAH .GT. 1> MRITE 16.12901 KSTART.KSTOP 

ha 

IN3T1 


CALL LOCATE IN. XPERT. YPERT. YCR2.10RAD.LMN2.LHX2.NCR2.LCR2.XLOC2) 

HAIN42S 



IF INPARAH .EO. It WRITE (6il300> 

HA 

IN372 


YMINaYPEBT(LMN2) 

MAIN429 



IF INPARAH .EO. It GO TO 210 

MA 

IN373 


YMAXaYPERT ILMX2) 

HAIN430 



NUHaKSTOP-KSTART*! 

HA 

1N37* 


YPCHINaYHiN 

HA1N431 



IF INUH .EO. 1) WRITE 16.13101 (OROIKI .KaKSTARTiKSTOPI 

ha 

IN375 


YPCHAXaYMAX 

HAIN432 



IF INUH .EO. 2) WRITE 16.1320) lOROIKI .KaKSTART.KSTOP) 

HA 

IN376 


WRITE 16.1380) ICASE.NCASE 

MAIN433 



IF INUH .EO. 3) WRITE 16.1330) lORDIK) .KaKSTART.KSTOP) 

MA 

IN37T 


WRITE 16.1090) VNAH 

HAIN434 



IF INUH .EO. 4) WRITE 16.13*0) lOROIK) .KaKSTART.KSTOP) 

HA 

IN370 


WRITE 16.1390) H2 

MAIN* 35 


210 

continue 

MA 

IN379 


IF INPARAH .EQ. 1) WRITE (6.1400) 02 1 1) .PARNAMI I ) 

MA1N436 



call fill (l.O.STRUNI) 

HA 

IN3P0 


IF INPARAH .GT. n WRITE 16.1410) (K.02 (K) .PARNAH(K) ,Ka| .NPARAH) 

HAIN437 



KLASTa204(KST0P*KSTART.H 

HA 

IN381 


WRITE 16.1130) VNAM.YCR2 

MAIN* 38 



WRITE I6.13S0) ISTRUNIlKI.Kal.KLAST) 

HA 

IN382 


WRITE (6.1140) 

MAIN* 39 



WRITE (6.1360) 

ha 

IN3K3 


WRITE (6.1150) 

HAIN4A0 



DO 220 lal.N 

HA 

1N3A* 


CALL UPLOW IA.B.XLOC2i6.NCR2*2.XOUT.FLAO) 

HAIN441 


220 

WRITE (6.1 370) I .XBASE 1 1 ) . IXUSAVE IK. | ) . YUN) T IK , I ) .KaKSTART .KSTOP ) 

HA 

IN38S 


WRITE 16.1420) HCAD2.X0UTIl).FLAGIl).LMN2.X0UTl2)«FL*G(2t.LHX2 

MAIN442 



IF lIRPT .EO. 0) GO TO 230 

HA 

IN386 


IF (NCR2 .GT. 01 WRITE (6.1430) NCR2. 

NA1N443 



IRPTaO 

HA 

1N38T 


« (OROIl) .X0UTI1*2).FLAG(I*2I .LCR2(I) .I«1.NCR2) 

MAIN444 



GO TO 190 

ma 

IN388 


CALL scale in. XBASE. 2.A.B) 

HAIN44S 


230 

call scale in. XBASE.). a. B) 

HA 

IN389 

c 


HA1N446 


2*0 

continue 

HA 

IN3O0 

c.... 

.IF LCHEX .NE. 0 READ IN DATA FOR COMPARISON SOLUTION ANO LOCATE 

Ha1N**7 

c 



HA 

IN391 

c 

MINIMUM. MAXIMUM. AND CRITICAL POINTS. 

HAtNAAa 

c. 

... 

CONSTRUCT PERTURBATION SOLUTIONS FOR TEST CASES (AND COMPARE WITH 

HA 

IN392 

c 


HAIN449 

c 


exact solution, if AVAILABLE). 

HA 

1N343 


IF ILCMEK .EO. 0) GO TO 290 

HA1N4S0 

c 



HA 

IN39* 


CALL INPUT IS) 

HA1N4S1 



00 330 ICASEal.NCASE 

MA 

1N39S 


call scale IN.xCHEK.I.A.B) 

HAIN4S2 



CALL INPUT U) 

ha 

INJ96 


CALL LOCATE (N,XCHeK.YCHEK,YCR3.IGH*D.LHN3.LMX3.NCR3.LCR3.XLOC3) 

HA1N453 



H2SAVEIlCASEIaM2 

HA 

1N397 


YMlNaYCHEKILMN3) 

HAIN4S4 



YCR2aYCRlTlH2) 

MA 

1N398 


YMAXaYCKEK (LMX3) 

MAINaSS 



YCR3aYC«2 

MA 



IF IY«IN .LT. yPCHIH) yPCHINxYHIM 

HAtSASfr 



IF (YM»X ,6T. VPCMAXI YPC»**X»VM*X 

CALL UPLOH tA»B.XL0C3»6tNCR3»?»X0UT.FLA0l 

WRlTt I6.U20> HEAO3tX0UT(I>.FLAGIIl.LMN3iXOUT(?I.FLAG12).LMX3 
IF (NCB3 .GT» 0> write (6»IA301 NCR3» 

% tOROU) »XOUT(I»?) .FLAGI I«?» tLCRS ( I ( . I » 1 »NCH3I 

CALL INTERP tN.XPERT.YPERT.XCHEK.YPRTIl 

CALL locate (NtXCH£K,YPRTI *YCR3« |GRAOtLMN3fLWX3«NCR3*LCR3tXLOC3) 
YTM1 N»yPRTI (LMN3) 

YTHAXxVPRTI ILHX 3I 
IF (YTMIN ,LT. YMIN) YMIN»YTMIN 
IF (YTHAX .GT. YHAX) YRAX«Y1MAX 
CALL SCALE (N.XPERT.J.AtBI 
call scale <NtXChEKt2«A«B) 

CALL FILL l2tO*STRING) 

WRITE (GalAAOl VNAW 

WRITE (GaUSO) VNAHfYMAXaVNAPaYMINaVNAMaTCRZaVNAMaVNAM 
WRITE iGaUGOt VNAMaWNAMaVNAMaVNABa (STRING! I ) a I*la72) 

00 200 I«laN 

call fill OalaSTRINGI 

200 WRITE IGaUTOI I aXBASE ( 1 ) a YBASE ( 1 1 aXPERT ( 1 1 a YPERT ( I ) a 
% XCHEK(I)aYCHEKCl) aYPRTI tI)a(STRING(IIlaII=laT2l 
GO TO 310 
290 CONTINUE 
C 

CALL SCALE INaXPERTa2aAaB> 

CALL FILL (2a0aSTRINGI 
WRITE (GalAAO) VNAN 

WRITE (GaUSO) VNAHaYHAXaVNANaYHINaVNAMaYCP2aVNAM 
WRITE (GalA90) VNAHa VNAMaSTRING 
DO 300 1>laN 
call FILL (3* I •STRING) 

300 WRITE (G.ISOO) I aXBASE ( 1 1 aYBASE 1 1 ) aXPERT ( 1 1 a YPERT ( I ) aSTRlNG 
310 continue 

IF LPLOT ,NE. 0 GENERATE PERIPHERAL PLOT OF PERTURBATION AND 

C COMPARISON solutions, 

C 

IF (LPLOT .EO. 0) 60 TO 320 

YMIN«Y0CHIN 

VMAXb VRCMAX 

IF (YPCHIN .LT. YMIN) YMIN«YPCH|N 
IF (TPCMAX .GT. YMAX) YNAX-YPCMAX 

CALL ORVPLT ( ICA5EaN,NCASEaNPARANaYNIN.YMAXaYCR2) 

320 continue 

CALL SCALE (NaXBASE a I a AaB) 

130 CONTINUE 
C 

PRINT FINAL TABLE S0NNAHI2IN6 CALCULATIONS. 

C 

call table (NPARAMaNCASEaPARNAMaHOaQOI 
GO TO 999 
C 

ABNORMAL TERMINATION OF COMPUTATION. 

C 

900 WRITE (Ga9000) 

GO TO 999 

90S WRITE (6a90S01 


MAlNASr 

MAINASa 

MAINAS9 

MAINAGO 

HAINAGl 

MAINAG2 

MAINAG3 

MAINAGA 

MAINAGS 

MAINAGG 

MAINAGT 

MAINAGS 

MAINAGO 

MAINA70 

MAINA71 

HAINA72 

MAINA73 

MAINA7A 

MAINA7S 

HAINA76 

MAINA77 

HAINA7B 

MAINA79 

MAINAAO 

HAINASI 

HAINAG2 

MAINA03 

MAINASA 

HAINASS 

HAINABG 

MAINA07 

MAINASS 

MAINAS9 

MAINA90 

MAINA9I 

MAINA92 

MAINA93 

MAINA9A 

HAINA9S 

MAINA96 

MAINA97 

MAINA90 

HAINA99 

MAINSOO 

MAINSOl 

MAINS02 

MAINS03 

HAINSOA 

MAINSOS 

MAINSOG 

HAINSO/ 

MAINSne 

MAINS09 

HAINSIO 

MAINSI I 

MAIN5I2 

MAINSI 3 


GO TO 999 
910 WRITE (6. 91001 
C 

999 WRITE (Ga95001 
STOP 
C 

I/O Format statements follow. 


1000 FORMAT 


1020 format 

s 

1030 format 


loso format 
lOGO format 
ioto format 
1080 format 
1090 format 
« 

% 

1100 FORMAT 
1110 format 
1120 format 
1130 format 
11 AO format 
USD format 
II GO format 
« 

1170 format 

« 

% 

1180 format 

9 

1190 FORMAT 
1200 FORMAT 
1210 FORMAT 
1220 FORMAT 

1210 format 

12A0 FORMAT 
12S0 FORMAT 
12G0 FORMAT 
1270 FORMAT 
1280 FORMAT 
« 

1S90 format 
1310 FORMAT 


(|HlaI32l|H«l/ 

lX,lH*a25X«80Ala2SXtlH*/ 

1X,132(IH»)///I 

llXa29H LIST OF INPUT PARAMETERS// 

GXaSHN «alA// 

GX.3HA >aFS.IaAXa3HB >aF5.|// 

GXaSHNPARAM «aI2///> 

(1X.22H STRAINING OPTIONS// 

GXa23MNUHBER OF FIXED POINTS! aI2/> 

IGXaAOHFIXED POINTS WILL BE PRESELECTED BY USER/ 
GXaSTHANO LISTED BELOW IN PRINTOUT FOR BASE/ 

GXa26HAN0 CALIBRATION SOLUTIONS.///) 

IGXaASHFIXEO POINTS WILL BE AUTOMATICALLY DETERMINED/ 
GXaAAHBV THE PROGRAM FOR ALL SOLUTIONS AS FOLLOWS)// 
lIKalAHTWO END POINTS) 

IllXalGHPOINT OF MINIMUM, IX,2A1 ) 

IllXalGHPOINT OF MAXIMUM. IX. 2AI ) 
lllXa2AlaGHCRlT < a AA.GHPOINT) I 
I1HI,26HRESULTS of COMPUTATIONS ON. IXaSAA///) 

lIXalTH MACH number./ 

6X,34HVALUES OF PERTURRATION PARAMETERS!/ 
GXalTHCRlTICAL VALUE OF a IXaZA] , IH I /) 

IlIXaAHHO ■aF7.A/> 

(llXaAHQO ■aF7.AaSXt|H(aA8a|H)/) 
lllXt3HOOI«I|«3H) >aFr.AaSXalH(.ABtlM)/) 

IllXf2AIa6HCRlT -.FBaA///) 

(1X.A7H LOCATIONS OF M|N.a MAX., ANO CRITICAL PTS.) 

13X,3AHI* denotes POINT ON LOWER SURFACE)/) 
IGXalAHHtNIMUM AT X ■aFT.AaAt a3Xa 1 OHIPOINT NO.alAalHI/ 
GX.IAHMAXIMUM AT X ■aFT.AaAl tax. lOHIPOlNT NO.alAalH)) 
(GX.II.IXalBHCRITICAL POINTIS)!/ 

(lOXaAAaGHAT X ■. 1X«F6. AaA| aJXa 

IGHIAFTER POINT NO.alAalH))) 

(lXaS(|H<) alXa32HFIXEO POINTS PRESELECTEO BY USER. 

IXaSIIH*)/) 

IGX,2HX(a|3alH)l 

(///1X.29H LOCATION OF FIXED POINTS) 

lllXa5HXFIX!aIla3H) ■aFT.AtAI) 

I1M1.2GHRESULTS OF COMPUTATIONS ON.lXaGAA///) 

IllXaAHMl >.F7.4/> 

l2Xa41H(«« DENOTES PERTURBATION FROM BASE VALUE)/) 
(I1X.4HOI ■aF7.4.5XalH(aABa|H)/> 
l9X,SH**01(alla3H) ■aFT.4.SX, IH I , ABa IH) /) 
!llXa3HQll«lla3H) «aFT.4,SX. 1H( ,A8a IH)/) 

!lH]a20HUNIT PERTURBATION OF . 1 X.2A I . IXa 
27HAND UNIT STRAINING OF XBASE) 

I2CBHF8R CALIBRATION SOLU1 IONS. I2.1X . THT)W)OUGH. 12) 
(///19Xal!|M«.lXaA4.|1HCALB SOLN w.SX)) 


MA 

HA 

NA 

MA 

MA 

MA 

MA 

HA 

HA 

MA 

HA 

HA 

MA 

MA 

HA 

MA 

HA 

MA 

HA 

MA 

MA 

MA 

MA 

HA 

NA 

MA 

MA 

HA 

HA 

HA 

MA 

HA 

MA 

MA 

MA 

HA 

MA 

MA 

MA 

HA 

MA 

MA 

ha 

NA 

HA 

MA 

HA 

MA 

HA 

HA 

MA 

MA 

MA 

MA 

MA 

HA 


INS|4 
INSIS 
INSIG 
INS IT 
INSIB 
1NS19 
INS20 
IN521 
IN522 
INS23 
INS24 
[NS2S 
INS28 
[NS27 
tNS28 
tN529 
IHS39 
IN531 
IN532 
INS33 
IN534 
IN53S 
INS3G 
IN537 
INS 38 
INS39 
INS40 
INSAl 
1NS42 
INS43 
[NS44 
INS4S 
[N54G 
IN54 7 
tNS48 
INS49 
[NSS8 
INSSl 
tN552 
tNSSJ 
INSS4 
[NSSS 
INSS6 
1N55T 
iNssa 
[N559 
[NSGO 
[NSGI 
INSG2 
INSG3 
INSG4 
[NSGS 
INSGG 
lNb67 
INS68 
|Nb70 



132* rOflHAT I///19X>2(1H*.1X«*4 i|1HCAL« SOLN •OX)) 

133fl FOWtXT l///I9X,3(|H«»lXtX*tUHCXLt SOLN •.3X)) 

1340 FONNAT l///|9Xi4(lH*<lX«*4»||HCALfl SOLN ••3X)) 

13S0 FORMAT (IXtSHPOlNT t4X«SHXBASE<4X >96X I ) 

1360 format MX) 

I3T0 FORMAT I |K, UtIXt I IF10.4) 

13X0 format llH1.32llH>t/ 

■ 1X.21H* OUTPUT FOR CASE N0.O3«4N OF tlZOH •/ 

• lXt32<lH*)///) 

1390 FORMAT lllXi4)W2 ««F7.4/) 

1400 FORMAT (1|Ai 4H02 *f FT.4»5X»1HI . A0> IHI /) 

1410 format MlX«3H02(ttl»3H) ••FTt4«SXi 1HUA6« IH) /) 

1420 format (/6Kt6A4// 

« llXiUHMINfMUM AT X itFT.AtAlOXilOHCPOINT N0.<I4«1H)/ 

« llXtl4HHAX|HUM AT X ■fF7.4iAl<3X< lOHIPOINT NU.fl4«IHI) 

1430 FORMAT <1H < 10X« tl« IXf 10)KRITICAL POINTIS)!/ 

« tl5XfA4»6HAT X ««F7.4iAl «3X« 

% I6HIAFTER POINT «0.<I4<1H) ) ) 

1440 FORMAT l///IX«44H... ..FINAL PRINTOUT AND GRAPHICAL DISPLAY OF.IX. 

S 2A|) 

1450 FORMAT </T2X.2lHN ■ MAXIMUM VALUE OF« |X »2A I . I IH> .FB. 4/ 

« T2X.21HL ■ minimum VALUE OF. IX . 2A 1 . I X . IH> .F S. */ 

» T2X.21H» ■ CRITICAL VALUE OF . 1* .2A I . IX. IH« ,FB . */ 

B T2X.12HP ■ value OF. I X.2A1 . IX. 

B 34HPREDICTC0 HY PERTURBATION SOLUTION/ 

B 72X.12HC > VALUE OF . 1X.2A1 . 1 X. 

B 22HIN COMPARISON SOLUTION/ 

B 72X.29HS 4 AGREEMENT BETWEEN P AND C) 

1460 FORMAT t/2X .2MPT .2X.5HXBASE.3X .2A1 .4HBA5E .2X.SHXPERT. 3X. 2A1 . 

B 4HPCRT.2X.SHXCHEK.3X.2A1.4HCHEX.2X.2A1.4HPINT, 1X.72AI/) 

1470 FORMAT ( IX . 1 3. 7F8.4. IX. 72AI > 

14B0 format </60X,21HH ■ MAXIMUM VALUE OF. IX.2A1 . t X « IM- .FS.4/ 

B 60X.21HL 4 minimum VALUE OF . | X >2A I . IX. |H« .F8.4/ 

B 60X«21H* 4 CRITICAL VALUE OF. |X .2A1 . IX. 1H< .FB.A/ 

B 60X.I2HP 4 VALUE OF . IX .2A1 . IX. 

B 34HPRED1CTC0 BY PERTURBATION SOLUTION) 

1490 format I/2X.2HPT.2X.5HXBASE.3X.2A1.4MBASC.2X.5HXPERT.3X.2A1. 

B 4HPERT.|X.96A|/) 

1500 FORMAT I IX. 1 3.4FB.4. 1X.96A1 ) 


abnormal TERMINATION FORMATS FOLLOW. 

C 

9000 format (///IX.2BHNUMSER OF CRITICAL POINTS IN/ 

B 1X.30HBASE AND CALIBRATION SOLUTIONS/ 

B 1X.3)MARE UNEQUAL - CALCULATION ENDED) 

9050 FORMAT (///IX.25HNUM8ER OF CRITICAL POINTS/ 

B IX.23HSELECTCn EXCEEDS NUMBER/ 

B 1X.3QHACTUALLY LOCATED - CALCULATION/ 

B IX.SHENOED) 

9100 FORMAT (///1X.2BHOROER OF SPECIFIED POINTS IN/ 

B lX,30MBASe AND CALIBRATION SOLUTIONS/ 

B 1X.39H00ES NOT CORRESPOND • CALCULATION ENDED) 

9500 format tlHD 
end 


SUBROUTINE BANNER 
C 

C PRINTS Two banner PAGES WITH JOB KEYWORD IN UPPER RIGHT CORNER. 

C 

COMMON /HEAD/ TITLE IBO) 

READ (S.IOOO) TITLE 
DO 10 tPAGC>1.2 

WRITE 16.1020) (TITLE(I).l4l.9) 

WRITE (6.1030) 

WRITE (6.1040) 

10 continue 

RETURN 

1000 format (80Al> 

1020 format (1H|.1I9X.I3(1H«)/120X«2H* .9A1.2H 4/120X. 13 ( ]H«) ) 

1030 format I////////38X.5S(1H*)/38X.1H».53X«1H*/ 

B 3XX.IH*. I9X.1SHPR0GRAM PERTURB. 19X . |H*/38X. |H*.S3Xi 1H»/ 

B 38X.IH*.7X.39HCALCULATES NONLINEAR MULTIPLE-PARAMETER. 

B TX.1H«/38X.1H».53X.1H«/ 

B 38X.1H4.13X.27HC0NTINU0US OR DISCONTINUOUS. 

B |3X.IH4/3BX.1M4,53x,1H»/ 

B 38X.1H*.1SX.22HPERTURBATI0N SOLUTIONS. 16X. 1H*/3BX. 1H*iS3X. 1H*/ 

B 38X.1H».10X.33HWHICH REPRESENT CHANGES IN EITHER. 

B |0X.1H*/3BX.1H<.S3X.1H*/ 

B 3BXi1H*.13X.27HGEONETRY OR FLOW CONDITIONS. 

B 13X.1H4/3BX.|H«.53X.1H*/ 

B 38X.IH*.4X.44HBY EMPLOYINO A STRAINEO'COOROINATE PROCEDURE. 

B 5X.1H«/3BX.IH»,53X.1H»/ 

B 38X.|H*i4X.44HUTILIZ1N0 UNIT PERTURBATIONS DETERMINED FROM. 

B 5X.IM»/38X.IH».53X.IH»/ 

B 3SX.|H«.I6X«2|HPREVI0USLY CALCULATED. 

B I6X.IH4/38X.1H4.53X.1H*) 

1040 FORMAT (3BX . IH4.9X.34K"BASER AND ’•CALIBRATION" SOLUTIONS. 

B I0X.1H»/38X.1H4.53X.1M»/ 

t 3eX.]H*.4X«45H0ISPLACC0 FRON ONE ANOTHER BT SOME REASONABLE. 

B 4X.IH»/3BX.1H4.S3X. |H*/ 

B 38X.|H>.SX.36HCHAN6E IN GEOMETRY OR FLOW CONDITION. 

B 9X.|M4/3bx,1H«,63X.1H4/ 

B 3BX. |H*.53Xi |H*/ 

B 3BX.IH4.2IX.10HWRITTEN BY .22X. 1H*/3BX« 1H>.53X. IH*/ 

B 38X.|H*i7x.39hJAMES P. ELLIOTT AND STEPHEN 5. STAHARA. 

B 7X.|H*/38X.1H*.S3X«|H*/ 

B 3nX.lH«.S3X.lH*/ 

B 38X.|H*.TX.38HNIELSEN ENGINEERING AND RESEARCH. INC.. 

B 8X. |H«/3BX,|h*.S3X.1H*/ 

B 3eX,lH«.UX,25HMOUHTAlN VIEW. CAL IFQHNl A, 14X» IH»/38X. IH».S3X. !«•/ 
B 3BX.55(1H»M 

END 


SUBROUTINE DRVPLT (ICA5C.N.NCA5C.NPARAM.YHIN.YMAX.VCR2) QRVP 1 

DIMENSION HLINE20) .HLINE30) ORVP 2 

DIMENSION XCSAVE(a.200).YCSAVE(8.200) ORVP 3 

DIMENSION XBASE<200).XCALB(200) .XPERT(?00) .XCHEKI200] «yBASE(200) . DHVP 4 

t YCALe<2001.YPERT«2081.YCHEK(200).YPRTt(2001 ORVP 5 




COMMON /HC«D/ ZTlTLEiaol 
COMMON /SAVE/ XCSAVEtYCSAVE 

COMMON /*y/ XO»SEt*CAL0tXHENT.*CHEK,Y0*5E»VC*LR.YPEBTfVCHtK.YPRTI 
DATA NPLOT /O/ 

IF (NPLOT ,EO. 01 CALL BETA 
MIN>in.0«<VM|N-0,l) 

MAX-10. 0*(YHAX*0.1) 

YMIN-0.I*MIN 

yhax-o.i*hax 

ENCODE (2S< lOOOtHLINEZ) (ZT I TLE ( t ) • t«l f 91 t ICASE .NCASL 

DO ZO K-1<NPARAM 

NPLOT.NPLOT-1 

ENCODE (22<1010tHLINE3> K.NPARAH 

CALL BGNPL 1-1 1 

CALL HIXALF ("L/CST") 

CALL MX3ALF ("INStR*S«%"l 
CALL SIMPLX 

CALL TITLE (IH . 1« IHXt I tM«E0.5l CtEX (P) S». I OOf 6.0.8. 0 1 
CALL HEADIN ("PLOT (OFl C*L0.25M0.7 (P) S". 100 .3.3) 

CALL HEADIN IHLINE2.28.2.3) 
call HEADIN IHLINE3.22.Z.3) 

call GRAF (0.0. "SCALE". I.O.VMAX. "SCALE". YMINI 

call frame 

IF (YCRZ .6T. VMAXI GO TO 10 
CALL RLVEC (O.O.YCRZ.0.2.YCR2.0000) 

CALL RLMESS ("C«L (P«BE (• I 1 00. 0 .21 . YCR2) 

10 CONTINUE 
CALL DASH 

CALL CURVE (XBASE. YBASE.N. 01 
CALL RESET ("DASH") 

call dot 

call copy (2.N.K.XCALB.XCSAVE) 
call copy (2.N.K.YCALB.YCSAVE) 

CALL CURVE (XCALB.YCALB.N.OI 
call reset ("DOT") 
call HANKER (1) 

CALL CURVE (XPERT. YPERT.N.-l) 

CALL RESET ("MARKER") 
call curve (XCHEK.YCHEK.N.O) 

CALL ENOPL (NPLOT) 

20 continue 

1000 FORMAT (9Al.)]Ht CASE NO. .12. AH OF .12) 

1010 FORMAT (16HCALIBRATI0N NO. .11. AH OF .11) 

RETURN 

end 


DRVP 6 
DRVP T 
DRVP B 
DRVP 9 
DRVP 10 
DRVP 11 
DRVP 12 
DRVP 13 
DRVP lA 
URVP IS 
DRVP 16 
DRVP IT 
DRVP 18 
DRVP 19 
DRVP 20 
DRVP 21 
DRVP 22 
DRVP 23 
DRVP 2A 
DRVP 2S 
DRVP 26 
DRVP 2T 
DRVP 28 
DRVP 29 
DRVP 30 
DRVP 31 
DRVP 32 
DRVP 33 
DRVP 3A 
DRVP 3S 
DRVP 36 
DRVP 3T 
DRVP 38 
DRVP 39 
DRVP aO 
DRVP Al 
DRVP a2 
DRVP A3 
DRVP AA 
DRVP AS 
DRVP A6 
DRVP AT 
DRVP a8 
DRVP 49 
DRVP SO 


SUflROUTINE ECHINP tCH| I 

dimension CARD(20) ECHl 2 

READ IS. 1000) LECHO ECHI 3 

IF (LECHO .EO. 0) GO TO 10 ECHl A 

WRITE (6.1010) ECHI 5 

BACKSPACE 5 ECHl 6 

BACKSPACE 5 ECHI T 

10 READ (S.I020) CARD ECHI B 

IF IC0F(5)| 30.20 ECHI 9 

2ft WRITE (1.1020) CARD ECHI |0 


IF (LECHO .EO. 1) WRITE (6.1030) CARO ECHI 11 

GO TO 10 ECHI 12 

30 REWIND 1 ECHI 13 

IF (LECHO .EO. 01 RETURN ECHl lA 

READ (1.1020) CARD ECHl IS 

READ (1.1020) CARD ECHl 16 

RETURN ECHI 17 

1000 format (151 ECHl 18 

1010 FORMAT (|H1.2S(1H»)/ ECHI 19 

* IX.1H*,1X.21HLISTING OF INPUT OECK.lX.IH-/ ECHI 20 

% 1X.2SIIH-)///) ECHI 21 

1020 FORMAT (20AA) ECHI 22 

1030 FORMAT (1X.20AA) ECHI 23 

end ECHI ?A 


SUBROUTINE FILL ( ICALL. I .STRING) FILL 1 

C f^lLL 2 

FILLS array STRING WITH CHARACTERS FOR TABLE HEADINGS AND PRINTER FILL 3 

C PLOTS. ^ILL a 

C I^ILL 5 

DIMENSION XBASEI200I.XCALB(200I.XCHEK(20a).XPERTI200) FILL 6 

DIMENSION YBASE(200).YCALB(200) .YPERT(200I .YCHEKI200) .YPRTI (200) FILL 7 

DIMENSION LOCO(6I .LOCI (6) .LSELCT (6) .00(8) .02(8) PILL 8 

DIMENSION STRINGI96).VNAHI2I.UN|TI20).PARNAMI8) FILL 9 

REAL MO. Ml. M2 FILL 10 

COMMON /PARAH/ PARNAH. LOCO. LOC I .LSELCT .N.NCASE .LSPEC.LUNI T . FILL II 

s lchekilplot.nselct.a.b.nparam.vnam fill 12 

COMMON /PERT/ 00.02.HO.M1.M2.OI.YCRO.YCR1.YCR2.YH1N.YMAX FILL 13 

COMMON /XY/ XBASE. XCALB. XPERT. XCHEK.YRASE.YCALR.YPERT.YCHEK.YPRTI FILL |A 

data I£NT /O/ pill 15 

DATA STAR/IM*/. P/IHP/. C/IHC/. OASH/IM-/. H/lHH/. EL/IHL/. PILL 16 

% 0LANX/IH /. DOLLAfl/lHS/ PILL 17 

DATA UNIT /IHX.IHS.IHT.IHR.IHU.IHN.IHI.IHT.IH .IH . FILL IB 

t IH .IH .IH .IHU.IHN.IHI.IHT.IH . IH , IH / FILL 19 

60 TO (10. 30. SO). ICALL PILL 20 

10 ient«ient*i pill 21 

IF (lENT .GT. II RETURN FILL 22 

UNIT(I2)>VNAH(1) PILL 23 

UNIT ( 131 «VNAH(2) FILL 2A 

DO 20 J-1.20 PILL 25 

DO 20 K-l.A PILL 26 

II»J*20*(K-1) PILL 27 

20 STP1NG(I1I»0NIT(J) PILL 28 

RETURN PILL 29 

3ft NY-96 PILL 30 

JFLAG-n PILL 31 

IF IYCR2 .GT. TNAX .OR. YCR2 .LT. YHIN) 1FLA6-I FILL 32 

IF (LCHEK .EO. 1) NT-72 PILL 33 

range-ymax-yhin pill 3A 

DO AO 11-1. NY PILL 35 

Aft STRING(1I)«DASH pill 36 

STRING(II-H PILL 37 

string (Ny) -EL PILL 38 

IF (IFLAG .EO. 1) RETURN PILL 39 

NSTAR>1*(VHAX-VCR2)/RAN6E*(NY-|I PILL AO 

string (NSIAR) -STAR PILL Al 

RETURN PILL ,A2 



SB continue fill *3 

DO 60 Il>l.NT fill 6* 

60 STRINGM1)>BL6NR FILL 45 

IF IIFLAG .EQ. 01 STRlNOINSTAR)«STAft FILL 46 

TP>Y(>ERT<1) fill 4T 

IF ILCHEK .EO. 1) YP-TPRTIII) FILL 48 

NPERT»|*(YM*X-YP1/R*NG£«(NY-1I fill 49 

STRING(NPERT)>P fill 50 

IF (LCHEK ,EO. 0) RETURN FILL 51 

YC»YCMEKIll FILL 52 

NCHEKs|*IYM6X-YC) /RANGE* INY-U FILL 53 

STRINGINCHEKKC fill 54 

IF (NPERT .EO. NCHEK) string (NPERT)sDOLLAR fill 55 

RFTUHN fill 56 

END FILL 5/ 


SUBROUTINE INPUT IICALLI INPU 1 
C INPU 2 
C WITH THE EXCEPTION OF THE TITLE AND THE ECHO PARAMETER 1 1 TENS I INPU 3 
C and 2)> ALL INPUT FOR PROGRAM PERTURB IS READ RT THIS SUBROUTINEt INPU 4 
C and is REOUIREO IN THE FOLLOMING ORDER IFOR DETAILS t REFER TO INPU 5 
C accompanying HANUALI. INPU 6 
r INPU T 
C***» ITEM NO. I - ONE CARO IBAIOI »••••••••••••••••••••••»•*••••••••••• INPU B 

C INPU 9 
C TITLE IDENTIFIES JOB - PRINTED AS HEADLINE ON FIRST PAGE INPU 10 
C OF OUTPUT. FIRST NINE CHARACTERS ARE PRINTED IN INPU || 
C UPPER RIGHT CORNER OF BANNER PAGE. AND IN UPPER INPU 12 
C LEFT CORNER OF SUMMARY PAGE. INPU 13 
C INPU |4 
C***» item no. 2 - ONE CARD II5> ••••••••**••••••••••••••••••••*•••••••• INPU 15 

r INPU 16 
C LECHO CONTROLS NHETHER OR NOT INPUT DECK IS PRINTED. INPU IT 
C INPU 18 
C LECHO ■ 0 ... NO OUTPUT INPU 19 
C LECHO ■ 1 ... OUTPUT INPU 20 
C INPU 21 
C***» ITEM NO. 3 - ONE CARD II6I5I ••••••*•*••••••••••••••*••••••••••••• INPU 22 

C INPU 23 
C N NUMBER OF DATA POINTS IN BASE AND CALIBRATION INPU 24 
C SOLUTIONS. INPU 25 
C INPU 26 
C NCASE NUMBER OF CASES FOR WHICH PERTURBATION SOLUTIONS ARE INPU 2T 
C TO BE COMPUTED. INPU 28 
C INPU 29 
C NPARAH NUMBER OF PARAMETERS PERTURBED. INPU 30 
C INPU 31 
C NSELCT number OF POINTS I IN ADDITION TO END POINTS) TO BE INPU 32 
C HELD INVARIANT IN STRAINING. INPU 33 
C NOTE) I .LE. NSELCT .LE. 6. INPU 34 
C INPU 35 
C LSPEC CONTROLS HOW INVARIANT POINTS IN STRAINING ARE INPU 36 
C SPECIFIED. INPU 3T 
C INPU 38 
C LSPEC ■ 0 ... INVARIANT POINTS SELECTED FROM AMONG INPU 39 
C those LOCATED BY THE PROGRAM INPu 40 
C <SEE ITEM NO. 4) INPU 4 1 


CO 


C 

C 

C 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

(•••• 

c 

c 

c 

c 
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c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c 

c**«* 

c 

c 

c 


LSPEC > 1 ... INVARIANT POINTS PRESELECTED BY USER INPU 42 
ISEE ITEMS NOS. 12 AND 16) INPU 43 

INPU 44 

LUNIT CONTROLS WHETHER OR NOT UNIT COORDINATE STRAININGIS) INPU 45 

AND UNIT PERTURBATION IS) ARC PRINTED. INPU 46 

INPU 47 

LUNIT ■ 0 ... NO OUTPUT INPU 48 

LUNIT * I ... OUTPUT INPU 49 

INPU 50 

LCHEK SPECIFIES WHETHER OR NOT PERTURBATION SOLUTION IS TO INPU 51 
BE CHECKED AGAINST AN EXACT COMPARISON SOLUTION. INPU 52 
A PRINTER PLOT OF SOLUTION IS MADE IN EITHER CASE. INPU 53 

INPU 54 

LCHEK ■ 0 ... NO COMPARISON INPU 55 

LCHEK • I ... COMPARISON INPU 5b 

INPU S7 

LPLOT SPECIFIES WHETHER OR NOT AN ADDITIONAL PLOT BY A INPU 58 

PERIPHERAL DEVICE IS TO BE MADE (SOFTWARE MUST BE INPU 59 

SUPPLIED BY USER IN SUBROUTINE PLOT). INPU 60 

INPU 61 

LPLOT > n ... NO PERIPHERAL PLOT INPU 62 

LPLOT > I ... PERIPHERAL PLOT INPU 63 

INPU 64 

item no. 4 - ONE CARD 1 1615) ••••••••••••••••••••••••*•••••••••••• INPU 65 

— this item to BE OMITTED IF LSPEC « 1 — INPU 67 

INPU 68 

LSELCTlI) ... INPU 69 

array of length 6 OF WHICH NSELCT ELEMENTS ARE READ INPU 79 

IN. SPECIFIES NATURE OF POINTS TO BE HELD INVARIANT INPU Tl 

ACCORDING TO THE CODE INPU T2 

INPU 73 

1 ... MINIMUM PT. HELD INVARIANT INPU 74 

2 ... maximum PT. held invariant INPU 75 

3 ... 1ST CRITICAL PT. HELD INVARIANT INPU 76 

4 ... END CRITICAL PT. HELD INVARIANT INPU 77 

5 ... 3RD CRITICAL PT. HELD INVARIANT INPU 78 

6 ... 4TH CRITICAL PT. HELD INVARIANT INPU 79 

INPU M8 

NOTE that the CODE NUMBERS CAN BE ASSIGNED IN ANY INPU 81 

ORDER. E.G. INPU 82 

INPU B3 

LSELCT(I) - 1 INPU M4 

LSELCTI2) - 3 INPU 85 

LSELCT(3> ■ 4 INPU 86 

INPU 87 

IS EQUIVALENT TO INPU 88 

INPU 89 

LSCLCT(l) • 4 INPU 99 

LSELCT(2) ■ 1 INPU 91 

LSELCTO) » 3 INPU 92 

INPU 93 

ROTH CORRESPONDING TO NSELCT « 3 WITH THE MINIMUM. INPU 94 

AND FIRST AND SECOND CRITICAL POINTS HELD INVARIANT. INPU 95 

INPU 96 

ITEM NO. 5 * ONE CARD I2AI) ••••••••••••••*•••••*•••*••••••••••»•• INPU 97 

INPU 98 

VNAM CHARACTER STRING OF LCN6TH 2 WHICH SVM80LIZCS INPU 99 

DEPENDENT VARIABLE. E.G. CP FOR PRESSURE 1NPU199 




o^ 


C COEFFICIENT. INPUIOl 

C INPUI02 

C*»*« ITEM MO. 6 - ONE CARO IIOARI *••*•••••••*••••••••••••••**•••*••••• IMPUI03 

C INPUIO* 

C PARNAM<KI ... INPUIOS 

C ARRAY OF 8-CHAMACTER STRINGS WHICH IDENTIFY THE INPUlOb 

C parameters varied. NPARAH ELEMENTS OF THE ARRAY INPUlOT 

C ARE READ IN. INPU108 

C INPUI09 

C»»*» item NO. 7 - ONE CAPO ISFlO.bl *••••»•••••••••••••••*••••••*••*••• INPUIIO 

C INPUllI 

C A scaling parameter (A « (I I . WHERE AID IS FIRST |NPU1|3 

C DATA POINT ON LOWER SURFACE ... SEE MANUAL). INPUII3 

C INPUIIA 

c S SCALINO PARAMETER IB > X|N>. WHERE AIN) IS LAST DATA INPUllS 

C POINT ON UPPER SURFACE ... SEE MANUAL). INPUlIb 

C INPUllT 

ITEM NO. • - ONE CARO CBFla.bl •••••••••••*••••••••••••••**••*«•** INPUU8 

C INPUU9 

C NO 0NC0MIN6 MACH NUMBER IN BASE SOLUTION. INPU120 

C INPU121 

item mo. 9 - ONE CARD I8FI0.6) INPU122 

C INPUI23 

C QOIK) array of length 8 OIVINO VALUES OF PERTURBATION INPUI2A 

C parameters in BASE SOLUTION. NPARAM ELEMENTS OF THE INPUI2S 

C ARRAY ARE READ IN. INPUI26 

C 1NPUI2T 

(•••• item no. 10 “ ONE SET OF C CAROS I8F10.6). C • I ♦ INTIN/8) *••••• 1NPU12S 

C INPUI29 

C XBASEIIlt 1>I.N ... INRUI30 

C A COORDINATE IN BASE SOLUTION. INPUI31 

C IMPUI32 

item NO. II - ONE set OF C CAROS I8FI0.6). C AS IN ITEM NO. 10 INPUI33 

C INPUI3A 

C VBASEII). lal.N ... INPUI39 

C OEPENOENT VARIABLE IN BASE SOLUTION. INPUI36 

C INPUI3T 

(•••• item no. 12 - ONE CARD I161S) 1NPU138 

C INPU139 

C this ITEM TO be omitted IF LSPEC ■ 0 INROUQ 

C INPUUl 

C LOCOII) ARRAY OF LENGTH 8 OF WHICH NSELCT ELEMENTS ARE READ INPUU2 

C IN. SPECIFIES SURSRIPTS OF THOSE BASE FLOW POINTS INPUU3 

c WHICH ARE TO BE HELD INVARIANT. INPUUA 

C 1NPUU5 

C*M«**»*»«**M«**»*»«»*»»*«*»»«***«««»»*«**»«**»**»»**«»»**»»*»*»«»*»»» INPUU* 
C INPUU7 

C ONE SET OF ITEMS 13 THROUGH 16 IS INPUUB 

C • NOTE * REQUIRED FOR EACH OF THE NPARAM INPUU9 

C •••••••• CALIBRATION SOLUTIONS. INPUlSO 

C INPUISI 

£•••• item NO. 13 - ONE CARO I8FI0.6I INPU152 

C INPU153 

C Ml ONCOMING MACH NUMBER IN KTH CALIBRATION SOLUTION. INPUISA 

C 1NPU1S5 

C 01 VALUE OF KTH PERTURBATION PARAMETER IN KTH INPUIS6 

C CALIBRATION SOLUTION. INPUIS7 


C INPUIS8 

C**«* ITEN NO. lA • ONE SET OF C CARDS ISF10.6)< C AS IN ITEM NO. 10 1NPU1S9 

C INPU160 

C XCALBIII. I* 1,N ... 1NPU16I 

C X COORDINATE IN KTH CALIBRATION SOLUTION. INPU162 

C 1NPUI63 

C»»** ITEM NO. IS - ONE SET OF C CAROS l8F10.6lt C AS IN ITEM NO. 10 INPUlbA 
C INPUlbS 

C YCALBIDt Iiil.N ... INPU166 

C dependent variable in KTH CALIBRATION SOLUTION,. INPUlbT 

r INPUlbS 

item no. lb - ONE CARO II61S) •••••*••••••*•••••»••••••••••••••••* INPlJlbR 

C INPUI70 

c THIS ITEM TO BE OMITTED IF LSPEC = 0 INPU171 

C IMPU1T2 

C LOCI (II ARRAY OF LENOTH 6 OF WHICH NSELCT ELEMENTS ARE READ INPUIT3 

C IN. SPECIFIES SUBSCRIPTS OF THOSE POINTS IN KTH INPUITA 

C CALIBRATION FLOW WHICH ARE TO BE MELD INVARIANT. INPUITS 

^ INPU176 

INPU17T 

C INPUIT8 

C ••«••••• ONE SET OF ITEMS |7 THROUGH 20 IS IMPU1T9 

C • note • REQUIRED FOR EACH OF THE NCASE INRUlBO 

C •••••••• SOLUTIONS TO BE COMPUTED. INPU18I 

C INPU182 

item no. IT - ONE CARD IBF10.6) •••••••••••••••••••••••••••*•••••• INPU183 

C INPUlBA 

C M2 ONCOMING MACH NUMBER IN SOLUTION TO BE COMPUTED. INPUISS 

C 1NRU1B6 

item no. 18 - ONE CARD I8F10.6) INPUI8T 

C INPUI88 

C 02IK1 ARRAY OF LENGTH B GIVING VALUES OF PERTURBATION INRUI89 

c parameters in SOLUTION TO BE COMPUTED. NPARAM INRU190 

C elements OF ARRAY ARE READ IN. INRU191 

C 1NPUI92 

ITEM NO. 19 - ONE SET OF C CAROS (8F 10.61. C AS IN ITEM NO. 10 ••• 1NPU193 

C INPU19A 

c — this ITEM TO BE OMITTED IF LCHEK ■ 0 INPUI95 

C INPUI96 

C XCHEKII), !■ I.N ... INPU19? 

C X COORDINATE IN COMPARISON SOLUTION. 1NPU19B 

C INPUI99 

item no. 20 - ONE SET OF C CARDS I8FI0.6). C AS IN ITEM NO. 10 ••• 1NPU200 

C INPU201 

C THIS ITEM TO BE OMITTED IF LCHEK > 0 INPU262 

C INPU203 

C VCHEKin, I>1.N ... INPU20A 

C DEPENDCNT VARIABLE IN COMPARISON SOLUTION. INPU20S 

C INPU206 

INPU207 

C INPU208 

DIMENSION LOCOI6I. LOCI 101. LSELCTI6I INPU209 

DIMENSION XBASEI200) IXCALBI200I iXPERT 1200) tXCHEK 12001 • INPU2I0 

% YBASE 1200) iVCALai200) .YPERT (200) .YCHEK 12001 .YPRTI <2001 INPU211 


DIMENSION 00 10) .0218) INPU212 
REAL MO.MI.MZ INPU213 
DIMENSION VNAMI2) .PARNAHIBI INPU2I4 



COMMON /P*R»M/ PARN*MfLOCO«LOCl«LSELCT,MfNCASE.LSPCC.LUNIT, 

% LCHCKiLP(.OTtNSELC1«*<BiNPAR«M.VN*M 

COMMON /PERT/ 0D*aZ>He>Hl.MZ.(ll,TCR0«<rCRl.YCR2.TMIN,TM«J( 

COMMON /XY/ X0»SEi*C*LB.KPERt.XCHEK»Yfl*SE.rCALB.YPEHT.YCHtKtYPRTI 
00 TO (l0«20>30t*D«50) • ICALL 

in READ Il»l000t H*MCASt.NPARAH,MSCLCT.l.SPEC»LUNJT»LCHEK.LPLOT 
IP CLSPCC .C0« 0) READ IlilOOO) ILSELCT ( 1 1 , 1 = 1 .NSELCT ) 

READ llttOlO) VNAM 

READ 11<10Z0> tPARNAM(K(«K>l<NPARAHI 

READ (WI030) AiB 

RETURN 

Zn READ ll>)»3«> MO 

READ <I«10301 lOO(K) >K*ltNPARAMI 
REAU l|<1030l (ABASE 1 1) f t = l (Nl 
READ Iltie30) lYBASEUI 

IF (LSPEC .EO. 11 READ (l«1000l (LOCO < 1 1 . 1*1 •NSC(-CT I 
RETURN 

30 READ (lfl030) HlfOl 

READ Ilil030> UCALBdl tl'IiN) 

READ Ut|030> (YCALe(I)»I*l»N> 

IF ILSPEC .ED. 1) READ UilOOOJ (LOCI ( ! I < I* 1 •MSELCT) 

RETURN 

40 READ <lil030t M2 

READ l!*1030) (02110 tK«l»HPAR AMI 
RETURN 

50 READ (|>I030) IXCHCK(fl.l>I«N) 

READ tl»103D) (YCMEMII 
RETURN 

1000 format IIOIB) 

IfllO FORMAT I2A1> 

1020 format (lOABI 
1030 format (BFIO.OI 
END 


INPU215 

INPU216 

1NPU21T 

INPU2I0 

INPU2I0 

INPU220 

INPU221 

INPU2Z2 

1NPU223 

INPU22A 

tNPU/?S 

INPUZZt. 

INRU/?7 

INPUZZH 

INPU229 

INPUZTO 

INPU23I 

INPU232 

INPU233 

INPU234 

INPU23S 

INPU236 

1NPU23T 

1NPU2JB 

INPU239 

INPU2A0 

INPU241 

INPU242 

1NPU2A3 

1NPU244 

1NPU2*S 

INPU240 

INPUZaT 


SUBROUTINE INTEBP (N«X,T.XI.Yn 

...GIVEN THE SET OF POINTS X(l»* Y(I)« I«1.N. 4N0 TmE SET XI (J). 

J>].Nt USES LINEAR INTERPOLATION TO COMPOTE THE SET YllJl* J»1»N. 

DIMENSION X(200l *Y (2001 .XH200I , YI (200» 

NM1>N-1 
JSTARTaI 
00 70 I*1«N 

IF (XKII .LE. Xdl) GO TO 10 
IF (XKlI .GE. X(Nlt GO TO 20 
GO TO 30 
in J»1 

GO TO GO 
20 J«M-1 
60 TO GO 
30 continue 

00 50 j*jstari,nmi 

IF (Xlll) .ME. X(J)) GO TO 40 

YI(1)>Y(J) 

GO rO 70 

40 IF (XKI) .GT. XUl .AND. XH1> .LT. X(J»1)I 60 TO GO 


INTE 

I 

INTE 

2 

INTE 

3 

INTE 

4 

INTE 

5 

INTE 

6 

INTE 

T 

INTE 

B 

INTE 

9 

INTE 

10 

JNTE 

11 

INTE 

12 

INTE 

13 

INTE 

14 

INTE 

IS 

INTE 

16 

INTE 

IT 

INTE 

IS 

INTE 

19 

INTE 

20 

INTE 

21 

INTE 

ZZ 


SO 

CONTINUE 

INTE 

23 

60 

SLOPE>(T(U*tl-V(JII/(X(U*ll-XIJM 

INTE 

24 


YI (IlnYljl •SLOPE* (XI <l)-KUII 

INTE 

?5 


JSTAHT«J 

INTE 

?6 

TO 

continue 

INTE 

27 


return 

INTE 

2S 


END 

INTE 

29 




subroutine LOCATE (NtX.ViYCRfr.rGRAD.LMrN.LM4X.NCRIT,LCRIT.XL0CI 

L0C4 

1 

C 



LOCA 

2 

c. 


.operates on the input array V. LOCATING MINIMUM AND MAXIMUM 

LOCA 

3 

c 


values* and all CRITICAL POINTS IY*TCR1T) FOR WHICH DY/OX (IN 

LOCA 

4 

c 


physical COOROINATESI has algebraic SION GIVEN BY IGNAD. NCH|T 

LOCA 

5 

c 


IS NUMBER OF CRITICAL POINTS. POINTS FOUND ARE STORED IN THE ARRAY 

LOCA 

6 

c 


XLOC AS FOLLOWS 

LOCA 

7 

c 



LOCA 

8 

c 


XLOCIll > MINIMUM PT. 

LOCA 

9 

c 


XL0CI2I * MAXIMUM PT. 

LOCA 

10 

c 


XL0CI3) • CRITICAL PT. NO. I 

LOCA 

11 

c 


... * • . . 

LOCA 

12 

c 


XL0C(6I > CRITICAL PT. NO. 4 

LOCA 

13 

c 



LOCA 

|4 



DIMENSION X (ZOO) .Y (200) <LCRIT (4) •XCRIT (4) «XLOC (6) 

LOCA 

IS 



COMMON /florev/ IREV 

LOCA 

16 



IFL0W*-1 

LOCA 

IT 



LMIN*1 

LOCA 

IB 



LH4X>I 

LOCA 

19 



ISTAHT*2 

LOCA 

20 



IF (IREV .EO. 01 GO TO 10 

LOCA 

21 



LMTN>z 

LOCA 

22 



LMAX-2 

LOCA 

23 



15TART»3 

LOCA 

Z4 


10 

CONTINUE 

LOCA 

25 



NCRIT>0 

LOCA 

26 



DO 30 I*ISTART.N 

LOCA 

2T 



IF (IREV .NE. 0 .AND. I .CO. N) GO TO 20 

LOCA 

2« 



IF (Y(II .GT. Y(LMAX>) LMAX*I 

LOCA 

29 



IF (Y(I) .LT. YILMIN)) LM1N*I 

LOCA 

30 


20 

CONTINUE 

LOCA 

31 



IF ((Y(I) ,GT. YCRIT .and. V(I«|I .61. YCRITi .OR. 

LOCA 

32 


( 

1 (Y(I1 .LT. YCRIT .AND. Y(l-l) .LT . YCRIT)) GO TO 30 

LOC* 

33 



IF (1 .GT. IREV) IFLOWiI 

LOCA 

34 



IF ( (Y(l)-Y(l-l))*FLOATItFLON*IGRAD) .LT. 0.0) GO TO 30 

LOCA 

35 



NCRITzNCRIT*! 

LOCA 

36 



LCRIT(NCRIT)*I-I 

LOCA 

37 



SL0PE«(XCI)-X(I-1))/(Y(I)-Y(I“I) ) 

LOCA 

3B 



XCRIT INCR1T)*X|1-1)»SL0PE*(YCRIT-YII-1>) 

LOCA 

39 


30 

continue 

LOC» 

*0 



XL0CII)*X(LM1N) 

LOCA 

41 



XLOC l?)«XILMAX) 

LOCA 

42 



IF (NCRIT .EO. 0) RETURN 

LOCA 

43 



DO 40 1*1.NCRIT 

LOCA 

44 


40 

XL0CII*Z)»XCH1T11) 

LOCA 

GS 



return 

LOCA 

46 



END 

LOCA 

47 


ID 



OO 




SUmtOUTtNC *WN0 IN.LtliTI 



MONO 

1 



SUBROUTINE SORT IN.X.ISEO) 

SORT 

I 

c 





MONO 

2 

C 



SORT 

2 

c« 


.ENFCK9 POINTS IN VICINITY Of A CRITICAL POINT FOR MONOTONE 


MONO 

3 

c. 


.ARRANGES THE SET XII). XI2I. ... . XIN) IN A MONOTONE INCREASING 

SOHT 

3 

c 


8EHAWtOR< AND ADJUSTS VALUES IF NECESSARY 

TO GIVE A LINEAR 


MONO 

4 

c 


SEQUENCE. ISEO GIVES ORDER OF SUBSCRIPTS IN REARRANGED SEQUENCE. 

SORT 

4 

c 


profile* 



MONO 

s 

c 



SOHT 

S 

c 





mono 

6 



DIMENSION XIOI.ISEQIS) 

SORT 

6 



DIMENSION LUI fXI200) •VI200I 



MONO 

7 



NMIaN-1 

SOHT 

7 



DO 10 I>1*N 



MONO 

a 



00 10 I>1.N 

SOHT 

a 



LS«un 



MONO 

9 


10 

ISFOII)>I 

SORT 

V 



Y1»Y<L5-U 



MONO 

10 


20 

ITESTaO 

SOHT 

10 



Y2«VILS) 



MONO 

11 



00 30 lal.NMl 

SOHT 

11 



V3aY<LS*I > 



MONO 

12 



IF ixin .le. xii*D) go to 3n 

SOHT 

12 



YA«T<LS*2> 



MONO 

13 



XSAVEaXlI) 

SOHT 

13 



IF MVI .LT* Y2I «AND. IY2 .LT. Y3> .AND. 

IY3 .LT. YAM GO 

TO in 

MONO 

lA 



XlllaXll.l) 

SOHT 

14 



IF IIVl .OT. V2I .ANO. <Y2 .GT. V3I .AND. 

IY3 .GT. YAM 60 

TO 10 

MONO 

IS 



XII.llaXSAvE 

SORT 

IS 



xisl rLS*H 



MONO 

10 



ISAVEalSEQII) 

SORT 

16 



X2aAILS> 



MONO 

IT 



ISEOCD-ISEOII*!) 

SORT 

17 



X3«X(LS*H 



MONO 

la 



ISE0II*H>ISAVE 

SORT 

IB 



XAal<LS*2> 



MONO 

IV 



ITESTal 

SORT 

IV 



SLOPE* «T4-T1»/«X*-X1) 



mono 

20 


30 

continue 

SOHT 

20 



T ILSI «Y 1 ‘SLOPE* < X2-X 1 1 



MONO 

21 



IF I1TE5T .EO. 1) 60 TO 20 

SORT 

21 



YILS»11«TI*SL0PE»«X3-X1I 



MONO 

22 



RETURN 

SORT 

22 


la 

CONTINUE 



MONO 

23 



ENO 

SORT 

23 



RETURN 



MONO 

?A 








ENO 



MONO 

2S 















SUBROUTINE strain IN.K.NSEO.XFIX.XIN.PaRM.OELX) 

stha 

I 








c 



STHA 

2 








c. 

. - .a 

•computes straining INCRENENT OELX from input array XIN. USING 

STRA 

3 



SUOROUTINE scale INiX«MtAtR> 



scal 

I 

c 


PIECEWISE LINEAR STRAINING WITH NSEG LINEAR SEGMENTS. FOR UNIT 

stha 

4 

c 





SCAL 

2 

c 


STRAINING, input VALUE OF PARM IS I.OI FOR GENERAL CASE. 

STHA 

s 

c. 


.ENTRY VITH M ■ 1 CONVERTS FROM PHYSICAL X 

10 TO -A ON LOWER 


SCAL 

3 

c 



stha 

6 

c 


SURFACE* 0 TO B ON UPPER SURFACE) TO N0RMALI2E0 A 10 .LT. X 

.LT. 

SCAL 

4 

c 


PARM a I02IK)-00(KI )/|OI>QOIK) ) . 

STHA 

7 

c 


It ENTRY VITH Ha? REVERSES THE PROCESS. NZ 

IOETERMINEO WHEN 

Hal) 

SCAL 

s 

c 



stha 

e 

c 


CORRESPONDS TO POINT AT NOSE OF BLADE OR AIRFOIL. 


scal 

6 



DIMENSION XFIX(8).XINI200) .OELX 1200) 

stha 

9 

c 





SCAL 

T 



COMMON /COEFF/ CtS.T) .OIS.T) 

STRA 

10 



COMMON /FLOREV/ NZ 



SCAL 

a 



JSTAHTal 

STHA 

11 



DIMENSION XI200) 



SCAL 

4 



DO 30 1*1. N 

STRA 

12 



IF (M .EO. 2> GO TO 30 



SCAL 

10 



00 10 JajSTART.NSEG 

STRA 

13 



CONTINUE 



SCAL 

11 



IF IXINIt) .6E. XFIXIJ) .AND. XINIII .LE. XFtXIJ‘1)) 60 TO 20 

STHA 

14 



NZ«0 



scal 

12 


10 

continue 

STMA 

IS 



00 10 I*2iN 



SCAL 

13 


20 

DELXII|aPARM*ICIK.JI‘(OIK.J)-I.O)*XlNIII ) 

STRA 

16 



IF (XU) .LT. XII'IM NZal 



scal 

1* 



JSTARTaJ 

STHA 

17 


10 

continue 



SCAL 

IS 


30 

CONTINUE 

STRA 

IB 



00 20 |al«N 



scal 

to 



RETURN 

STHA 

IV 



IF (I .LE. NZ) Ta-Xll) 



SCAL 

17 



end 

STRA 

20 



IF It ,OT. NZI TaXdl 



scal 

)B 








XII)a|T-.At/(0-A) 



SCAL 

IV 







20 

CONTINUE 



SCAL 

20 



SUBROUTINE TABLE INPARAN.NCASE.PARNAM.MO.QOI 

TABL 

1 



RETURN 



SCAL 

21 



DIMENSION PARNAMIBI.QOIB) 

TAOL 

2 


30 

00 *0 lal.N 



SCAL 

22 



dimension 0RDI2S) .HEADOIS) .HEAD) IS) .HEA02ISI .MISAVE IB) • 

tabl 

3 



X(I)aAeS(<B-AI*X(II*A) 



SCAL 

23 


S M2SAVEI2S) .OISAVEI0) .02SAVE IR.2S) 

TABL 

4 


40 

continue 



SCAL 

24 



dimension owriteioi 

TABL 

5 



RETURN 



SCAL 

25 



REAL H0.N1SAVE.M2SAVE 

lAHL 

6 



end 



SCAL 

26 



COMNON /HEAD/ TITLE 100) 

TABL 

7 



COMHON /TW./ HC*D*tHC«Dl«HC«D2.N|S«VE«H2S*VE*0IS«VEt02SAVE TAtfL 8 

OlHENStON STAR (1271 fBLANK 1 1261 (SrPNeK 100) tSTRMC2l28)iENOI6) 7*6L V 

DATA star /127«|H</* SCANK /126*1H /, STRN02 / 1M« ,27» IH-/ . TAHL 10 

t STRMOl /100*1H-/ TABL 11 

DATA END /IHPf IHAtlHRtlR.tlH tlH*/ TABL 12 

DATA ORO /5H 1ST tSH 2ND .5H 3RD .5H ATM «SH STM tSM 6TM tSH 7TM f TABL 13 

% 5M BTH fSH »TH tSMlOTH tSHUTM .SH12TH .SH13TM .SMIATH i TABL lA 

% 5H1STH tSH16TH tSMlTTH tSHIBTH .SH19TH •SM20TH i5H2|5T » TABL IS 

S 5H22ND •SH23RD fSH2ATM .5H2STM / TABL 16 

00 10 I»1.10 TAHL IT 

II«11»I-10 TABL 18 

10 STRNG1(II»«STARI11 TAHL 19 

00 50 IPAGE«1»2 TABL 20 

WRITE I6tl000l ITITLEII>»I-1«9I TABL 21 

1HAA1>11*NRARAH*12 TABL 22 

IMAA2>IHAX|*28 TABL 23 

WRITE I6tl0101 IBLANKtI>tI-l«28)«ISTAR<I)fl>l*IMAXI) TABL 2A 

WRITE I6tl020> lORDIKItKvltNPARAMl TABL 25 

IHAX>IMAX2>6 TAHL 26 

WRITE l*»103ll (BLANK ( II flaltlNAXI •END TABL ?T 

IMAXaIMAXl-11 TABL 28 

WRITE (6«10A0I ISTRN«ltI><taI<|HAXI TABL 29 

WRITE ( 6 « 10501 (PARNAM(K)«RaltNPARAN| TABL 30 

]MAKaIHAX2-I TABL 31 

WRITE (6tl0301 (M.ANKIII«Ial>IHAXIfSTAR(H TABL 32 

WRITE (6*1010l (STARIIIfIaI*IMAX2l TABL 33 

WRITE (6tlOAOI HEADOtNOt(M(KltK>I.NPARAHI TaBL 3A 

WRITE <6«1030l (BLANK(IliI>ltIHAXItSTAR(l) TaBL IS 

WRITE (6tl010> (STAR(|lfI>If|HAX2l TABL 36 

WRITE (6«1070) (STAR(II«I«1>1MAX2) TABL 37 

DO 30 KaltNPARAH TABL 38 

DO 20 KKaltNPARAN TABL 39 

to OWRlTE(K)OaOO(KKI TABL AD 

OWRITECKIaOlSAVEIKI TABL Al 

WRITE (6.10801 ORO(K).HEADttHlSAVE(KI.(QWRITEIKK),KKa|.NPARAHI TABL A2 

WRITE (6.10301 (BLANK I n.I-l.IMAX). STAR! 1) TABL A3 

IF (K .LT. NPARAHt WRITE 16.10101 STRNG2. (STRNGl ( 1 1 . la ] . IHAXl I TABL AA 


30 continue TABL as 

WRITE (6.10101 (STARdl .lal.lMAXt) TABL A6 

WRITE (6.10701 (STARdl .|al.|MAX2) TABL a7 

DO AO ICASE-I.NCASE TABL a8 

WRITE (6.10801 0RD(ICASEI.MEA02.M2SAVEdCASE) . TABL A9 

% IQ2SAVE(K.ICASE).Kal,l4PARAMI TABL SO 


WRITE (6.10301 (BLANKIII.lal. IRAK). STARCH TABL SI 

IF dCASE .LT. NCASEI WRITE (6.10101 STRNG2. (STRNGl ( 1 1 . 1 = 1 . IMAX I ) TABL S2 

AO continue TABL S3 

WHITE (A.lOlOl (STARdl .Ia].IMAX2l TABL 5A 

SO CONTINUE TABL 5S 

RETURN TABL S6 


1000 FORMAT ( IMl . AX . 1 31 IN*! /SX.2H* .9A1.2M */SX. I3C |M*) ///] 
1010 FORMAT <5X.12BA1I 

1020 FORMAT l33X.lMa.|0I.BdH9.A5.5MPAR. II 
1030 FORMAT I1H..4X.I28A1I 
lOAO FORMAT I33X.I1H* MACH NO. .89AH 
lOSO FORMAT l 33 X.lHa.ltX. 0 ( 2 M* .AB.IXll 


TABL ST 
TABL 58 
TABL S9 
TABL 60 
TABL 61 
TABL 62 


to 

to 


I 


1060 FORMAT (SX.3H* .SAA.SX.OOh* .F6.3.2XI1 TABL 63 

1070 FOPHAT (1H0.AX.I2BAII TABL 6A 

1080 format (SX.2M* . A5.5AA. 1 X.9( 3M* .F6.3.2XII TABL 6S 

ENO TABL 66 



subroutine UPLOW (A.8.XIN.K.N.X0UT.FLAGI 

UPLO 

1 

c 


UPLO 

2 

c. 

....CONVERTS NORMALIZED ARRAY XIN TO PMTSICAL ARRAY XOUT AND FLAGS 

UPLO 

3 

c 

POINTS ON LOWER SURFACE WITH A •••". 

UPLO 

A 

c 


UPLO 

5 


OIHENSION XINIKI.X0UTI8I 

UPLO 

6 


DIMENSION FLAGI8I 

UPLO 

7 


DATA BLANK/IM /. STAR/IH*/ 

UPLO 

8 


XNOSEa-A/IB'AI 

UPLO 

9 


DO 10 lal.N 

UPLO 

10 


FLAGdIaBLANK 

UPLO 

11 


IF (XINdl .LT. XNOSEI FLAGIIIaSTAR 

UPLO 

12 


XOUTIIIaABSdB-AiaXlNdl .Al 

UPLO 

13 


10 CONTINUE 

UPLO 

lA 


RETURN 

UPLO 

IS 


end 

UPLO 

16 


I 



APPENDIX C 
LIST OF SYMBOLS 



blade chord, m 

design variable coefficient of profile shape 
function; eq. (22) 

invariant point index; eq. (5); also, index for 
surface shape functions; eqs . (22,23) 

dummy index; eq . (24) 

two-dimensional full potential operator; eq. (1) 

linear operator representing first-order perturbation 
of two-dimensional full potential equation; eq. (3) 

linear operator representing first-order perturbation 
terms arising from coordinate straining; eq. (6) 

number of independent flow or geometrical variables 
to be perturbed 

absolute inlet Mach number 


total number of shock points and high-gradient 
maxima points; eq. (21) 

total number of invariant points, equal to n + 2 ; 
eqs . (13,14) 

til 

j arbitrary geometric or flow parameter to be 
perturbed; eq. (8) 


calibration flow value of q^ ; eq. (8) 
base flow value of q^. ; eq. (2) 

approximate flow solution for arbitrary flow quantity; 
eq. (7) 


calibration flow solution for value q of arbitrary 
parameter; eq . (7) j 


base flow solution for values q of arbitrary 
parameters; eq. (7) j 

t h 

j perturbation solution per unit change of 
perturbed parameter q^. ; eq. (7) 


100 


Cs,t) 

t 

Cx,y) 






Strained (x,y) coordinates; eq. (4) 
gap to chord spacing ratio 

nondimens ional blade-fixed orthogonal coordinates; 
eq. C7) , normalized by C 

nondimens ional blade-fixed orthogonal coordinates 
related to j calibration solution; eq. (7) 

straining functions associated with (x,y) coordinates; 
eq. (4) 


Cx-j ,y-| ) straining functions associated with ith invariant 
■^i "^i point; eq. (5) 

CL angle oncoming flow makes with blade chord line 


C6x^,6y^) 


e . 


J 


e 


j 


unit displacements in (x,y) directions associated 
with ith invariant point; eqs . (5,8) 

t h 

desired perturbation change of j geometric or 
flow parameter; eq. (8) 

perturbation change of j geometric or flow parameter 
between base and calibration flows; eq. (8) 


T 


thickness ratio of blade 


$ 


4 


o 


4 


nondimensional total velocity potential; eq. (1), 
normalized by CV 

nondimensional base flow velocity potential; eq . (2), 

normalized by CV 

t h 

nondimensional j perturbation velocity potential; 
eq. (2), normalized by CV^ 


Subscripts 

i denotes quantities associated with i invariant 

point 

j denotes perturbation quantities 


Superscripts 

o denotes base flow quantities 

c denotes quantities associated with calibration flows 
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TABLE 1 


COMPARISON OF FINAL DESIGN VARIABLES AND OBJECTIVE FUNCTION 
WHEN EMPLOYING FULL NONLINEAR TSONIC SOLUTIONS OR 
PERTURBATION METHOD FOR DIFFERENT CHOICES OF 
CALIBRATION SOLUTION MATRIX FOR SIX DESIGN 
VARIABLE SUBCRITICAL OPTIMIZATION CASE 
STUDY USING MAXIMUM SUCTION SURFACE 
VELOCITY DIFFUSION OBJECTIVE 


Design 


Objective 

Variables 

KOCR 

1 T 

ZM 

1 P 

TMX 

THLE 

Function 

INITIAL 

Baseline 
Upper Bound 
Lower Bound 

-10.0000 

0.0000 

-15.0000 

0.2500 

0.6000 

0.2000 

0.4500 

0.5500 

0.2000 

1.50000 

4.0000 

0.5000 

0.0500 

0.1000 

0.0300 

0.0050 

0.0120 

0.0030 

1.8400 




FINAL 







TSONIC 

SOLUTIONS 

ONLY RESULTS 



CDC 7600 
IBM 3033 

-7.1106 

-8.8659 

0.2000 

0.2400 

0.5500 

0.5500 

0.9401 

0.7628 

0.0300 
0.0300 _ 

0.0064 

0.0051 

1.6748 

1.6752 

PERTURBATION SOLUTION RESULTS 

CASE 1 

Calibration 

Final 

-7.0000 

-9.2223 

0.2000 

0.2327 

0.5500 

0.5500 

0.9400 

0.9281 

0.0300 

0.0359 

0.0064 

0.0052 

1.6904 

CASE 2 

Calibration 

Final 

-9.0000 

-8.8714 

0.2300 

0.2228 

0.5500 

0.5500 

0.8000 

0.9523 

0.0300 

0.0315 

0.0060 

0.0051 

1.6908 

CASE 3 

Calibration 

Final 

-12.0000 

-8.9865 

0.3000 

0.3001 

0.4000 

0.5500 

1.2500 

0.8776 

0.0400 

0.0300 

0.0070 

0.0050 

1.6974 

CASE 4 

Calibration 

Final 

-8.0000 

-8.2180 

0.3500 

0.3257 

0.5000 

0.5500 

2.0000 

1.4527 

0.0600 

0.0385 

0.0060 

0.0049 

1.7628 

CASE 5 

Calibration 

Final 

-9.0000 

-5.9055 

0.3000 

0.4440 

0.5000 

0.5500 

2.5000 

1.3325 

0.0400 

0.0412 

0.0040 

0.0037 

1.7486 

CASE 6 

Calibration 

Final 

-11.0000 

-9.4297 

0.2300 

0.2522 

0.3500 

0.5500 

1.2500 

0.8036 

0.0600 

0.0300 

0.0040 

0.0052 

1.6807 
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Figure 1.- Comparison of perturbation (O) and nonlinear ( — ) 
surface pressures for the simultaneous two-parameter 
perturbation of ^or nonlifting strongly 

supercritical flows past isolated 
NACA OOXX blade profiles. 
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OBJECTIVE IS TO MINIMIZE; 



+ 0.5 



Figure 5.- Illustration of physical basis of optimization problem involving 
blade surface contouring to taylor the surface pressure distribution 

to a desired distribution. 




o Full Aero, 

• Pert. Method 
■ Final Aero. Check 



Figure 6.- Comparison of final design variables and objective function when 
employing the perturbation method (•) in lieu of the full nonlinear 
aerodynamic solver (O) for various choices of initial design 
variable stepsize for nine design variable subcritical 
optimization case study with a surface pressure 
tailoring objective. 






Figure 7.- Comparison o£ computational work and objective 
function reduction per optimization search cycle when 
employing perturbation method (•) or full nonlinear 
aerodynamic solver nine design variable 

subcritical optimization case study using a 
surface pressure tailoring objective. 
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Figure 8.- Comparison o£ perturbation-predicted final design variables 
objective function for various choices of initial design variable 
stepsize for four design variable supercritical optimization 
case study with drag minimization objective. 


and 






Figure 9.- Comparison o£ computational work and objective function 
reduction per optimization search cycle when employing 
perturbation method (•) or full nonlinear aerodynamic 
solver Co) for four design variable supercritical 
optimization case study using a drag 
minimization objective. 
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Figure 10.- Comparison of final design variables and objective function when 
employing perturbation method (•) or full aerodynamic solver (A) as 
compared with 'optimum' full aerodynamic result (0) for 
different choices of initial design variable stepsize 
for four design variable supercritical optimization 
case study using a surface pressure 
tailoring objective. 



Figure 11.- Comparison o£ computational work and objective 
reduction per optimization search cycle when employing 
perturbation method (•) or full nonlinear 
aerodynamic solver (©) for four design 
variable supercritical optimization 
case study using a surface pressure 
tailoring objection. 
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